diff --git a/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/Advancing_front_surface_reconstruction.txt b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/Advancing_front_surface_reconstruction.txt new file mode 100644 index 00000000000..cee26494920 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/Advancing_front_surface_reconstruction.txt @@ -0,0 +1,36 @@ + +namespace CGAL { +/*! + +\mainpage User Manual +\anchor Chapter_Advancing_Front_Surface_Reconstruction +\anchor I1ChapterAdvancingFrontSurfaceReconstruction +\cgalAutoToc +\author Tran Kai Frank Da and David Cohen-Steiner + +\image html afsr.png + +Assume we are given a set \f$ S\f$ of points in 3D and we'd like to + + +\section AFSR_Definitions Definitions + +We distinguish + +\section AFSR_Examples Examples + +\subsection AFSR_Example_function Example for Global Function + +The basic +\cgalExample{Advancing_front_surface_reconstruction/reconstruction_fct.cpp} + + +\subsection AFSR_Example_class Example for Class + +A + +\cgalExample{Advancing_front_surface_reconstruction/reconstruction_class.cpp} + +*/ +} /* namespace CGAL */ + diff --git a/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction.h b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction.h new file mode 100644 index 00000000000..e06b4c26dd0 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction.h @@ -0,0 +1,152 @@ + +namespace CGAL { + +/*! +\ingroup PkgAdvancingFrontSurfaceReconstruction + +The class `Advancing_front_surface_reconstruction` + +\tparam Dt must be a `Delaunay_triangulation_3` with +`Advancing_front_surface_reconstruction_vertex_base_3` and `Advancing_front_surface_reconstruction_cell_base_3` +blended into the vertex and cell type. + + +\cgalHeading{Implementation} + +The .. +*/ +template< typename Dt> +class Advancing_front_surface_reconstruction { +public: + +/// \name Types +/// @{ + +/*! + The type of the 2D triangulation data structure describing the reconstructed surface. + The type `TDS_2::Vertex` is model of the concept `TriangulationDataStructure_2::Vertex` and has additionally the + method `vertex_3()` that returns a `Triangulation_3::Vertex_handle` to the associated 3D vertex + The type `TDS_2::Face` is model of the concept `TriangulationDataStructure_2::Face` and has additionally the + method `facet()` that returns the associated `Triangulation_3::Facet`, and a method `bool is_on_surface()` + for testing if a face is part of the reconstructed surface or a face incident to a boundary edge. + In case the surface has boundaries, the 2D surface has one vertex which is associated to the infinite + vertex of the 3D triangulation. +*/ + typedef Hidden_type TDS_2; + +/*! +The type of the 3D triangulation. + +*/ + typedef Dt Triangulation_3; + + +/*! + A bidirectional iterator which allows to enumerate all points that were removed + from the 3D Delaunay triangulation during the surface reconstruction. The value type + of the iterator is `Triangulation_3::Point_3`. +*/ +typedef Hidden_type Outlier_iterator; + +/*! + A forward iterator which allows to visit all boundaries. It + visits the entry point of each boundary twice. This allows to + detect that the traversal of a boundary is finished. One more increment + brings us to the vertex on the next boundary. + The value type of the iterator is `Triangulation_3::Vertex_handle`. +*/ +typedef Hidden_type Boundary_iterator; + + +/// @} + +/// \name Creation +/// @{ + +/*! +Initializes from a 3D Delaunay triangulation of a point set. +*/ +Advancing_front_surface_reconstruction(Dt& dt); + + +/// @} + +/// \name Operations +/// @{ + +/*! +calls the surface reconstruction function with the default parameters. +*/ + void operator()(); + +/*! +returns the reconstructed surface. +*/ +const TDS_2& +tds_2(); + +/*! +returns the underlying 3D Delaunay triangulation. +*/ +const Triangulation_3_& +triangulation_3(); + + + + +/*! +An iterator over the outliers. +*/ +Outlier_iterator outliers_begin(); + +/*! +Past-the-end iterator. +*/ +Outlier_iterator outliers_end(); + +/*! +An iterator over the boundary vertices. +*/ +Boundary_iterator boundary_begin(); + +/*! +Past-the-end iterator. +*/ +Boundary_iterator boundary_end(); + +/// @} + +/// \name Predicates +/// @{ + +/*! +returns `true` if the reconstructed surface has boundaries. +*/ +bool +has_boundaries() const; + + +/*! +returns `true` if the facet is on the surface. +*/ +bool +has_on_surface(Triangulation_3::Facet f) const; + +/*! +returns `true` if the facet f is on the surface. +*/ +bool +has_on_surface(TDS_2::Face_handle f2) const; + +/*! +returns `true` if the facet f is on the surface. +*/ +bool +has_on_surface(TDS_2::Vertex_handle v2) const; +/// @} + + + + +}; /* end Advancing_front_surface_reconstruction */ +} /* end namespace CGAL */ diff --git a/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction_cell_base_3.h b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction_cell_base_3.h new file mode 100644 index 00000000000..3b902c4e79d --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction_cell_base_3.h @@ -0,0 +1,22 @@ + +namespace CGAL { + +/*! +\ingroup PkgAdvancingFrontSurfaceReconstruction + +The class `Advancing_front_surface_reconstruction_cell_base_3` is the default +cell type for the class `Advancing_front_surface_reconstruction`. + +\tparam Traits has to be a model of `DelaunayTriangulationTraits_3`. + +\tparam Cb has to be a model of `TriangulationCellBase_3`. + +*/ +template< typename Traits, typename Cb = Triangulation_cell_base_3 > +class Advancing_front_surface_reconstruction_cell_base_3 : public Cb { +public: + +/// @} + +}; /* end Advancing_front_surface_reconstruction_cell_base_3 */ +} /* end namespace CGAL */ diff --git a/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction_vertex_base_3.h b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction_vertex_base_3.h new file mode 100644 index 00000000000..49414059b97 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/CGAL/Advancing_front_surface_reconstruction_vertex_base_3.h @@ -0,0 +1,23 @@ + +namespace CGAL { + +/*! +\ingroup PkgAdvancingFrontSurfaceReconstruction + +The class `Advancing_front_surface_reconstruction_vertex_base_3` is the default +vertex type for the class `Advancing_front_surface_reconstruction`. + +\tparam Traits has to be a model of `DelaunayTriangulationTraits_3`. + +\tparam Vb has to be a model of `TriangulationVertexBase_3`. + + +*/ + template< typename Traits, typename Vb = Triangulation_vertex_base_3 > +class Advancing_front_surface_reconstruction_vertex_base_3 : public Vb { +public: + +/// @} + +}; /* end Advancing_front_surface_reconstruction_vertex_base_3 */ +} /* end namespace CGAL */ diff --git a/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/PackageDescription.txt b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/PackageDescription.txt new file mode 100644 index 00000000000..2f53214e304 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/PackageDescription.txt @@ -0,0 +1,27 @@ +/// \defgroup PkgAdvancingFrontSurfaceReconstruction Advancing Front Surface Reconstruction Reference +/// \defgroup PkgAdvancingFrontSurfaceReconstructionConcepts Concepts +/// \ingroup PkgAdvancingFrontSurfaceReconstruction + +/*! +\addtogroup PkgAdvancingFrontSurfaceReconstruction + +\cgalPkgDescriptionBegin{Advancing Front Surface Reconstruction,PkgAdvancingFrontSurfaceReconstructionSummary} +\cgalPkgPicture{afsr-detail.png} +\cgalPkgSummaryBegin +\cgalPkgAuthors{Tran Kai Frank Da, David Cohen-Steiner} +\cgalPkgDesc{This package offers a surface reconstruction algorithm.} +\cgalPkgManuals{Chapter_Advancing_Front_Surface_Reconstruction,PkgAdvancingFrontSurfaceReconstruction} +\cgalPkgSummaryEnd +\cgalPkgShortInfoBegin +\cgalPkgSince{2.1} +\cgalPkgDependsOn{\ref PkgTriangulation3Summary} +\cgalPkgBib{cgal:dc-afsr} +\cgalPkgLicense{\ref licensesGPL "GPL"} +\cgalPkgShortInfoEnd +\cgalPkgDescriptionEnd + +This chapter presents a ... The description is based on +the articles \cite toto. + +*/ + diff --git a/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/examples.txt b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/examples.txt new file mode 100644 index 00000000000..2f14a4e739c --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc/Advancing_front_surface_reconstruction/examples.txt @@ -0,0 +1,4 @@ +/*! +\example Advancing_front_surface_reconstruction/reconstruction_fct.cpp +\example Advancing_front_surface_reconstruction/reconstruction_class.cpp +*/ diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction/main.aux b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction/main.aux new file mode 100755 index 00000000000..01003caf766 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction/main.aux @@ -0,0 +1,34 @@ +\relax +\@writefile{toc}{\contentsline {chapter}{\numberline {1}Advancing Front Surface Reconstuction}{1}} +\@writefile{lof}{\addvspace {10\p@ }} +\@writefile{lot}{\addvspace {10\p@ }} +\@writefile{lof}{\contentsline {xchapter}{Advancing Front Surface Reconstuction}{1}} +\@writefile{lot}{\contentsline {xchapter}{Advancing Front Surface Reconstuction}{1}} +\newlabel{chap:surface_reconstruction}{{1}{1}} +\@writefile{toc}{\contentsline {section}{\numberline {1.1}Introduction}{1}} +\@setckpt{Advancing_front_surface_reconstruction/main}{ +\setcounter{page}{2} +\setcounter{equation}{0} +\setcounter{enumi}{0} +\setcounter{enumii}{0} +\setcounter{enumiii}{0} +\setcounter{enumiv}{0} +\setcounter{footnote}{0} +\setcounter{mpfootnote}{0} +\setcounter{part}{0} +\setcounter{chapter}{1} +\setcounter{section}{1} +\setcounter{subsection}{0} +\setcounter{subsubsection}{0} +\setcounter{paragraph}{0} +\setcounter{subparagraph}{0} +\setcounter{figure}{0} +\setcounter{table}{0} +\setcounter{r@tfl@t}{0} +\setcounter{LT@tables}{0} +\setcounter{LT@chunks}{0} +\setcounter{mtc}{1} +\setcounter{minitocdepth}{2} +\setcounter{ptc}{0} +\setcounter{parttocdepth}{2} +} diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction/main.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction/main.tex new file mode 100755 index 00000000000..713ddc1e4fe --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction/main.tex @@ -0,0 +1,13 @@ +\chapter{Advancing Front Surface Reconstuction} +\label{chap:surface_reconstruction} +\ccChapterAuthor{Frank Da and David Cohen-Steiner} + + +\section{Introduction} + +This package offers an algorithm for surface reconstruction from an unorganized point set. +The algorithm selects facets of the 3D Delaunay triangulation of the points. + + + + diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_cell_base_3.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_cell_base_3.tex new file mode 100755 index 00000000000..d7cc8686c18 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_cell_base_3.tex @@ -0,0 +1,37 @@ +\begin{ccRefClass}{AFSR_cell_base_3} + +\ccDefinition + The class \ccRefName\ is the cell class that must be blended in the cell parameter +of the 3D Delaunay triagulation. It mainly provides storage and methods used by the +Advancing front surface reconstruction algorithm. + + +\ccInclude{CGAL/AFSR_cell_base_3.h} + +\ccInheritsFrom + +\ccc{CellBase} + +%\ccTypes + +%\ccTypedef{typedef Kernel::FT FT;}{The number type.} + +%\ccCreation +%\ccCreationVariable{c} + +%\ccConstructor{AFSR_cell_base_3();}{Default constructor} + +%\ccOperations + +%\ccMemberFunction{bool has_facet_on_surface(int i) const;}{Returns \ccc{true}, iff the facet is on the surface.} + +%\ccSeeAlso + +%\ccRefIdfierPage{} + +\end{ccRefClass} + + + + + diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_options.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_options.tex new file mode 100755 index 00000000000..06b2ffae373 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_options.tex @@ -0,0 +1,14 @@ +\begin{ccRefClass}{AFSR_options} + +\ccDefinition + +The class \ccRefName\ is used as a container of the options of the surface reconstruction algorithm. + +\ccInclude{CGAL/AFSR_options.h} + + +%\ccSeeAlso + +%\ccRefIdfierPage{} + +\end{ccRefClass} diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_3.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_3.tex new file mode 100755 index 00000000000..bffba82de7d --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_3.tex @@ -0,0 +1,33 @@ +\begin{ccRefClass}{AFSR_vertex_base_3} + +\ccDefinition + +The class \ccRefName\ is the vertex class that must be blended in the vertex parameter +of the 3D Delaunay triagulation. It mainly provides storage and methods used by the +Advancing front surface reconstruction algorithm. + +\ccInclude{CGAL/AFSR_vertex_base_3.h} + + +\ccInheritsFrom + +\ccc{VertexBase} + +%\ccTypes + +%\ccTypedef{typedef Kernel::FT FT;}{The number type.} + +%\ccCreation +%\ccCreationVariable{v} + +%\ccConstructor{Surface_vertex_base_2();}{Default constructor} + +%\ccOperations + +%\ccMemberFunction{Vertex_handle_3 vertex_3() const;}{Returns the vertex handle in the 3D Delaunay triangulation.} + +%\ccSeeAlso + +%\ccRefIdfierPage{} + +\end{ccRefClass} diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_with_id_3.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_with_id_3.tex new file mode 100755 index 00000000000..276a2aef984 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_with_id_3.tex @@ -0,0 +1,31 @@ +\begin{ccRefClass}{AFSR_vertex_base_with_id_3} + +\ccDefinition + +The class \ccRefName\ is the vertex class that must be blended in the vertex parameter +of the 3D Delaunay triagulation. It mainly provides storage and methods used by the +Advancing front surface reconstruction algorithm. + +\ccInclude{CGAL/AFSR_vertex_base_with_id_3.h} + + + +%\ccTypes + +%\ccTypedef{typedef Kernel::FT FT;}{The number type.} + +%\ccCreation +\ccCreationVariable{v} + +%\ccConstructor{Surface_vertex_base_2();}{Default constructor} + +\ccOperations + +\ccMemberFunction{int& id() const;}{Returns a reference to an \ccc{int} that is not used or altered by the + surface reconstruction algorithm.} + +%\ccSeeAlso + +%\ccRefIdfierPage{} + +\end{ccRefClass} diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Advancing_front_surface_reconstruction.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Advancing_front_surface_reconstruction.tex new file mode 100755 index 00000000000..6f405227aab --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Advancing_front_surface_reconstruction.tex @@ -0,0 +1,71 @@ +\begin{ccRefClass}{Advancing_front_surface_reconstruction} + +\ccDefinition + +The class \ccRefName\ extracts a surface from a 3D Delaunay triangulation. + +\ccInclude{CGAL/Advancing_front_surface_reconstruction.h} + +\ccParameters + +The parameter for \ccc{Delaunay_3} must be a 3D Delaunay triangulation where \ccc{CGAL::AFSR_vertex_base_3} +or \ccc{CGAL::AFSR_vertex_base_with_id_3} and \ccc{CGAL::AFSR_cell_base_3} must be blended in the vertex and face class. + + + +\ccTypes + +\ccTypedef{typedef Delaunay_3 Triangulation_3;}{} +\ccTypedef{typedef CGAL::Triple Edge;}{For \ccc{(ch,i,j)}, this is the edge between vertices \ccc{i} and \ccc{j} in cell \ccc{*ch}.} + +\ccTypedef{typedef std::pair Edge_incident_facet;}{For \ccc{((ch,i,j),k)}, this is the facet adjacent to the edge \ccc{(,i,j)}, and + opposite to vertex \ccc{k}, in the cell \ccc{*ch}.} +\ccGlue +\ccNestedType{TDS_2}{The type of the 2D triangulation data structure describing the reconstructed surface.} +\ccNestedType{TDS_2::Vertex}{It is model of the concept \ccc{TriangulationDataStructure_2::Vertex} and has additionally the + method \ccc{vertex_3()} that returns a \ccc{Triangulation_3::Vertex_handle} to the associated 3D vertex.} +\ccNestedType{TDS_2::Face}{It is model of the concept \ccc{TriangulationDataStructure_2::Face} and has additionally the + method \ccc{facet()} that returns the associated \ccc{Triangulation_3::Facet}.} + +\ccNestedType{Boundary_iterator}{This forward iterator allows to visit all contours. It + visits the entry point of the contour twice. This allows to + detect that the traversal of the border is finished. One more increment + brings us to the vertex on the next boundary. The value type of the iterator is \ccc{Triangulation_3::Vertex_handle}.} + +\ccNestedType{Outlier_iterator}{This bidirectional iterator allows to enumerate all points that were removed + from the 3D Delaunay triangulation during the surface reconstruction. The value type + of the iterator is \ccc{Kernel::Point_3}.} + +\ccCreation +\ccCreationVariable{es} + +\ccConstructor{Advancing_front_surface_reconstruction(Delaunay_3& del, AFSR_options opt);}{} + +\ccOperations + +\ccMemberFunction{const TDS_2& tds_2() const;}{Returns a const reference to the reconstructed surface. } +\ccMemberFunction{bool is_on_surface(TDS_2::Vertex_handle vh;}{Returns true, iff the vertex is on the reconstructed surface.} + + +\ccMemberFunction{const Triangulation_3& triangulation() const;}{Returns a const reference to the 3D Delaunay triangulation.} + +\ccMemberFunction{Boundary_iterator boundaries_begin() const;}{} +\ccMemberFunction{Boundary_iterator boundaries_end() const;}{} + +\ccMemberFunction{Outlier_iterator outliers_begin() const;}{} +\ccMemberFunction{Outlier_iterator outliers_end() const;}{} + +\ccMemberFunction{bool has_on_surface(Triangulation_3::Facet f) const;} + {Returns \ccc{true}, iff the facet is on the reconstructed surface.} + +\ccMemberFunction{Edge_incident_facet next(const Edge_incident_facet& f) const;}{returns the next facet around the edge.} +\ccMemberFunction{Edge_incident_facet previous(const Edge_incident_facet& f) const;}{returns the previous facet around the edge.} + +\ccMemberFunction{Facet next_surface_facet(const Edge_incident_facet& f) const}{returns the next facet ariund the edge, + which is on the reconstructed surface.} + +%\ccSeeAlso + +%\ccRefIdfierPage{} + +\end{ccRefClass} diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Surface_face_base_2.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Surface_face_base_2.tex new file mode 100755 index 00000000000..d9328eba201 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Surface_face_base_2.tex @@ -0,0 +1,39 @@ +\begin{ccRefClass}{AFSR::Surface_face_base_2} + +\ccDefinition + +The class \ccRefName\ is the face class used in the triangulation data structure +that describes the reconstructed surface. A face is either part of this +surface, or it is not part of the surface but adjacent to the boundary of the surface. + +The face class stores a facet of the 3D Delaunay triangulation, if it +is part of the surface. + +\ccInclude{CGAL/AFSR/Surface_face_base_2.h} + + +%\ccTypes + +%\ccTypedef{typedef Kernel::FT FT;}{The number type.} + +%\ccCreation +\ccCreationVariable{f} + +%\ccConstructor{Surface_face_base_2();}{Default constructor} + +\ccOperations + +\ccMemberFunction{Facet facet() const;}{Returns the facet in the 3D Delaunay triangulation.} + +\ccMemberFunction{bool is_on_surface() const;}{Returns \ccc{true}, iff the face is on the surface.} + +%\ccSeeAlso + +%\ccRefIdfierPage{} + +\end{ccRefClass} + + + + + diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Surface_vertex_base_2.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Surface_vertex_base_2.tex new file mode 100755 index 00000000000..2a2845e231d --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/Surface_vertex_base_2.tex @@ -0,0 +1,29 @@ +\begin{ccRefClass}{AFSR::Surface_vertex_base_2} + +\ccDefinition + +The class \ccRefName\ is the vertex class used in the 2D triangulation data structure +that describes the reconstructed surface. It stores a vertex handle to its +correpsonding vertex in the 3D Delaunay triangulation. + +\ccInclude{CGAL/AFSR/Surface_vertex_base_2.h} + + +%\ccTypes + +%\ccTypedef{typedef Kernel::FT FT;}{The number type.} + +%\ccCreation +\ccCreationVariable{v} + +%\ccConstructor{Surface_vertex_base_2();}{Default constructor} + +\ccOperations + +\ccMemberFunction{Vertex_handle_3 vertex_3() const;}{Returns the vertex handle in the 3D Delaunay triangulation.} + +%\ccSeeAlso + +%\ccRefIdfierPage{} + +\end{ccRefClass} diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/intro.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/intro.tex new file mode 100755 index 00000000000..ee184dd930d --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/intro.tex @@ -0,0 +1,16 @@ +\chapter{Advancing Front Surface Reconstruction} +\label{chap:surface_reconstruction_ref} +\ccChapterAuthor{Frank Da and David Cohen-Steiner} + + +\ccHeading{Classes} + +\ccRefIdfierPage{CGAL::Advancing_front_surface_reconstruction}\\ +\ccRefIdfierPage{CGAL::AFSR_options}\\ +\ccRefIdfierPage{CGAL::AFSR_face_base_3}\\ +\ccRefIdfierPage{CGAL::AFSR_vertex_base_3}\\ +\ccRefIdfierPage{CGAL::AFSR_vertex_base_with_id_3} + +\ccHeading{Functions} + +\ccRefIdfierPage{CGAL::write_to_file_vrml2} diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/main.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/main.tex new file mode 100755 index 00000000000..1501aa839d0 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/main.tex @@ -0,0 +1,8 @@ + +\input{Advancing_front_surface_reconstruction_ref/intro.tex} + +\input{Advancing_front_surface_reconstruction_ref/Advancing_front_surface_reconstruction.tex} +\input{Advancing_front_surface_reconstruction_ref/AFSR_cell_base_3.tex} +\input{Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_3.tex} +\input{Advancing_front_surface_reconstruction_ref/AFSR_vertex_base_with_id_3.tex} +\input{Advancing_front_surface_reconstruction_ref/write_to_file_vrml2.tex} \ No newline at end of file diff --git a/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/write_to_file_vrml2.tex b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/write_to_file_vrml2.tex new file mode 100755 index 00000000000..cc7ffd4a287 --- /dev/null +++ b/Advancing_front_surface_reconstruction/doc_tex/Advancing_front_surface_reconstruction_ref/write_to_file_vrml2.tex @@ -0,0 +1,18 @@ +\begin{ccRefFunction}{write_to_file_vrml2} + +\ccDefinition + +The function \ccRefName\ writes the reconstructed surface to a file. + +\ccInclude{CGAL/IO/Advancing_front_surface_reconstruction.h} + +\ccFunction{template +void +write_to_file_vrml2(char* foutput, const Surface& S, + bool boundary, double red, double green, double blue, bool no_header);} +{Opens a file with the basename \ccc{foutput} and the file extension \ccc{.wrl}, and writes the +reconstructed surface to this file. If \ccc{boundary} is \ccc{true} the boundary is drawn. +The color of the facets is defined by \ccc{red}, \ccc{green}, and \ccc{blue}. If \ccc{no_header} is +\ccc{true}, the VRML header is not written.} + +\end{ccRefFunction} \ No newline at end of file diff --git a/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/CMakeLists.txt b/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/CMakeLists.txt new file mode 100644 index 00000000000..1ff756f5b94 --- /dev/null +++ b/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/CMakeLists.txt @@ -0,0 +1,33 @@ +# Created by the script cgal_create_cmake_script +# This is the CMake script for compiling a CGAL application. + + +project( Advancing_front_surface_reconstruction_example ) + +CMAKE_MINIMUM_REQUIRED(VERSION 2.4.5) + +set(CMAKE_ALLOW_LOOSE_LOOP_CONSTRUCTS true) + +if ( COMMAND cmake_policy ) + cmake_policy( SET CMP0003 NEW ) +endif() + +find_package(CGAL QUIET COMPONENTS Core ) + +if ( CGAL_FOUND ) + + include( ${CGAL_USE_FILE} ) + + include( CGAL_CreateSingleSourceCGALProgram ) + + include_directories (BEFORE ../../include) +# include_directories (BEFORE ../../../../../trunk/Triangulation_2/include) + + create_single_source_cgal_program( "extract.cpp" ) + +else() + + message(STATUS "This program requires the CGAL library, and will not be compiled.") + +endif() + diff --git a/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/extract.cpp b/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/extract.cpp new file mode 100644 index 00000000000..5543d0d0e41 --- /dev/null +++ b/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/extract.cpp @@ -0,0 +1,468 @@ +#define NOLAZY +#define BLIND + + +#include +#include +#include +#include + +#include + +// Kernel +#include +#include +#include + +#include +#include +#include +#include +#include +#include + + + + + +typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; + +typedef Kernel::Point_3 Point; +typedef Kernel::Vector_3 Vector; + +typedef CGAL::AFSR_vertex_base_with_id_3 LVb; + +typedef CGAL::Triangulation_cell_base_3 Cb; +typedef CGAL::AFSR_cell_base_3 LCb; + +typedef CGAL::Triangulation_data_structure_3 Tds; +typedef CGAL::Delaunay_triangulation_3 Triangulation_3; + +typedef Triangulation_3::Vertex_handle Vertex_handle; + +typedef CGAL::Advancing_front_surface_reconstruction Surface; +typedef CGAL::AFSR_options Options; + + +//--------------------------------------------------------------------- + +struct Auto_count : public std::unary_function >{ + mutable int i; + Auto_count() : i(0){} + std::pair operator()(const Point& p) const { + return std::make_pair(p,i++); + } +}; + + +bool +file_input(const Options& opt, std::vector& points) +{ + const char* finput = opt.finname; + bool xyz = opt.xyz; + + std::ios::openmode mode = (opt.binary) ? std::ios::binary : std::ios::in; + std::ifstream is(finput, mode); + + if(opt.binary){ + CGAL::set_binary_mode(is); + } + if(is.fail()) + { + std::cerr << "+++unable to open file for input" << std::endl; + exit(0); + return false; + } + else + std::cerr << "Input from file : " << finput << std::endl; + + std::size_t n; + if(! xyz){ + is >> n; + std::cerr << " reading " << n << " points" << std::endl; + points.reserve(n); + CGAL::cpp11::copy_n(std::istream_iterator(is), n, std::back_inserter(points)); + } else { + // we do not know beforehand how many points we will read + std::istream_iterator it(is), eof; + char ignore[256]; + while(it!= eof){ + points.push_back(*it); + is.getline(ignore,256); // ignore what comes after 3 doubles in a line + it++; + } + n = points.size(); + } + + return true; +} + + + +void usage(char* program) +{ + std::cerr << std::endl << "NAME " << std::endl + << program << " - surface extension -" << std::endl << std::endl; + + std::cerr << std::endl << "OPTIONS" << std::endl + << " -xyz : input data in xyz format" << std::endl + << " -no_border -nb : " << std::endl + << " -in fname : reads points from file ./fname" << std::endl + << " -out fname : writes points to file ./fname" << std::endl + << " -out_format -of : choose file format for output (iv, wrl, off, medit," << std::endl + << " ply, stl, all, none)" << std::endl + << " -rgb r g b : color of the surface" << std::endl + << " -no_header : The Vrml header and footer are not written" << std::endl + << " -area a : No faces larger than area * average_area" << std::endl + << " -perimeter p : No faces larger than perimeter * average_perimeter" << std::endl + << " -abs_area a : No faces larger than abs_area" << std::endl + << " -abs_perimeter p : No faces with perimeter longer than abs_perimeter" << std::endl + << "\n Options for internal use" << std::endl + + << " -sect_in fname : reads points from sections file ./fname" << std::endl + << " -binary : binary I/O" << std::endl + << " -delta x : set the delta constant" << std::endl + << " -ki x y : set the K interval (default : [1.1 5])" << std::endl + << " -ks x : set the K step (default : .1)" << std::endl + << " -k x : set the K constant (only one pass)" << std::endl + << " -Delaunay : display the underlying Delaunay triangulation" << std::endl + << " -max_of_connected_components x : set the max of connected components" << std::endl + << " (default : non-active)" << std::endl + << " -post x : set a number for the post process" << std::endl + << " -contours : display contours" << std::endl; +} + + + +bool +parse(int argc, char* argv[], Options &opt) +{ + std::strcpy(opt.program, argv[0]); + --argc; + argv++; + if(argc == 0) + std::cerr << "nothing ???" << std::endl; + + while ((argc > 0) && (argv[0][0] == '-')){ + if ((!std::strcmp(argv[0], "-D")) || (!std::strcmp(argv[0], "-Delaunay"))) { + opt.Delaunay = true; + argv++; + argc--; + std::cerr << "-D "; + } + else if ((!std::strcmp(argv[0], "-c")) || (!std::strcmp(argv[0], "-contours"))) { + opt.contour = true; + argv++; + argc--; + std::cerr << "-c "; + } + else if ((!std::strcmp(argv[0], "-b")) || (!std::strcmp(argv[0], "-binary"))) { + opt.binary = true; + argv++; + argc--; + std::cerr << "-b "; + } + else if ((!std::strcmp(argv[0], "-x")) || (!std::strcmp(argv[0], "-xyz"))) { + opt.xyz = true; + argv++; + argc--; + std::cerr << "-x "; + } + else if ((!std::strcmp(argv[0], "-nb")) || (!std::strcmp(argv[0], "-no_border"))) { + opt.K = HUGE_VAL; + opt.K_init = opt.K; + argv++; + argc--; + std::cerr << "-nb "; + } + else if ((!std::strcmp(argv[0], "-nh")) || (!std::strcmp(argv[0], "-no_header"))) { + opt.no_header = true; + argv++; + argc--; + std::cerr << "-nh "; + } + else if ((!std::strcmp(argv[0], "-d")) || (!std::strcmp(argv[0], "-delta"))){ + if (sscanf(argv[1], "%lf", &opt.delta) != 1) { + std::cerr << "Argument for delta must be a number" + << std::endl; + } + argv += 2; + argc -= 2; + std::cerr << "-d " << opt.delta << " "; + } + else if ((!std::strcmp(argv[0], "-a")) || (!std::strcmp(argv[0], "-area"))){ + if (sscanf(argv[1], "%lf", &opt.area) != 1) { + std::cerr << "Argument for area must be a number" + << std::endl; + } + argv += 2; + argc -= 2; + std::cerr << "-a " << opt.area << " "; + } + else if ((!std::strcmp(argv[0], "-pe")) || (!std::strcmp(argv[0], "-perimeter"))){ + if (sscanf(argv[1], "%lf", &opt.perimeter) != 1) { + std::cerr << "Argument for perimeter must be a number" + << std::endl; + } + argv += 2; + argc -= 2; + std::cerr << "-perimeter " << opt.perimeter << " "; + } + else if ((!std::strcmp(argv[0], "-aa")) || (!std::strcmp(argv[0], "-abs_area"))){ + if (sscanf(argv[1], "%lf", &opt.abs_area) != 1) { + std::cerr << "Argument for abs_area must be a number" + << std::endl; + } + argv += 2; + argc -= 2; + std::cerr << "-abs_area " << opt.abs_area << " "; + } + else if ((!std::strcmp(argv[0], "-ae")) || (!std::strcmp(argv[0], "-abs_perimeter"))){ + if (sscanf(argv[1], "%lf", &opt.abs_perimeter) != 1) { + std::cerr << "Argument for abs_perimeter must be a number" + << std::endl; + } + argv += 2; + argc -= 2; + std::cerr << "-abs_perimeter " << opt.abs_perimeter << " "; + } + else if ((!std::strcmp(argv[0], "-ki"))){ + if ((sscanf(argv[1], "%lf", &opt.K_init) != 1)|| + (sscanf(argv[2], "%lf", &opt.K) != 1)){ + std::cerr << "Argument for K must be a number" + << std::endl; + } + argv += 3; + argc -= 3; + std::cerr << "-ki " << opt.K_init << " " << opt.K << " "; + } + else if ((!std::strcmp(argv[0], "-rgb"))){ + if ((sscanf(argv[1], "%lf", &opt.red) != 1)|| + (sscanf(argv[2], "%lf", &opt.green) != 1) || + (sscanf(argv[3], "%lf", &opt.blue) != 1)){ + std::cerr << "Argument for rgb must be three numbers" + << std::endl; + } + argv += 4; + argc -= 4; + std::cerr << "-rgb " << opt.red << " " << opt.green << " " << opt.blue << " " ; + } + else if ((!std::strcmp(argv[0], "-ks"))){ + if (sscanf(argv[1], "%lf", &opt.K_step) != 1) { + std::cerr << "Argument for K must be a number" + << std::endl; + } + argv += 2; + argc -= 2; + std::cerr << "-ks " << opt.K_step << " "; + } + else if ((!std::strcmp(argv[0], "-k"))){ + if (sscanf(argv[1], "%lf", &opt.K) != 1) { + std::cerr << "Argument for K must be a number" + << std::endl; + } + opt.K_init = opt.K; + argv += 2; + argc -= 2; + std::cerr << "-k " << opt.K_init << " "; + } + else if ((!std::strcmp(argv[0], "-m")) || (!std::strcmp(argv[0], "-max_of_connected_components"))){ + if (sscanf(argv[1], "%d", &opt.max_connected_comp) != 1) { + std::cerr << "Argument for the number of connected components must be a number" + << std::endl; + } + /* + if(opt.max_connected_comp < 1) { + std::cerr << "Argument for the number of connected components must be a positive number" + << "It is set to 1" << std::endl; + opt.max_connected_comp = 1; + } + */ + argv += 2; + argc -= 2; + std::cerr << "-m " << opt.max_connected_comp << " "; + } + else if ((!std::strcmp(argv[0], "-p")) || (!std::strcmp(argv[0], "-post"))){ + if (sscanf(argv[1], "%d", &opt.NB_BORDER_MAX) != 1) { + std::cerr << "Argument for post process must be a number" + << std::endl; + } + argv += 2; + argc -= 2; + std::cerr << "-p " << opt.NB_BORDER_MAX << " "; + } + else if ((!std::strcmp(argv[0], "-i")) || (!std::strcmp(argv[0], "-in"))) { + std::strcpy(opt.finname, argv[1]); + opt.file_input = true; + argv += 2; + argc -= 2; + std::cerr << "-i " << opt.finname << " "; + } + else if ((!std::strcmp(argv[0], "-s")) || (!std::strcmp(argv[0], "-sect_in"))) { + std::strcpy(opt.finname, argv[1]); + opt.Section_file = true; + opt.file_input = true; + argv += 2; + argc -= 2; + std::cerr << "-s " << opt.finname << " "; + } + else if ((!std::strcmp(argv[0], "-o")) || (!std::strcmp(argv[0], "-out"))) { + std::strcpy(opt.foutname, argv[1]); + opt.file_output = true; + argv += 2; + argc -= 2; + std::cerr << "-o " << opt.foutname << " "; + } + else if ((!std::strcmp(argv[0], "-of")) || (!std::strcmp(argv[0], "-out_format"))) { + if (!std::strcmp(argv[1], "wrl")) + opt.out_format = 0; + else if (!std::strcmp(argv[1], "off")) + opt.out_format = 1; + else if (!std::strcmp(argv[1], "medit")) + opt.out_format = 2; + else if (!std::strcmp(argv[1], "ply")) + opt.out_format = 3; + else if(!std::strcmp(argv[1], "iv")) + opt.out_format = 4; + else if(!std::strcmp(argv[1], "stl")) + opt.out_format = 5; + else if (!std::strcmp(argv[1], "all")) + opt.out_format = -1; + else if (!std::strcmp(argv[1], "none")) + opt.out_format = -2; + else + std::cerr << "unrecognized file format." << std::endl; + opt.file_output = true; + std::cerr << "-of " << argv[1] << " "; + argv += 2; + argc -= 2; + } + else if ((!std::strcmp(argv[0], "-?")) || + (!std::strcmp(argv[0], "-h")) || + (!std::strcmp(argv[0], "-help"))) { + usage(opt.program); + return false; + } + else { + std::cerr << "unrecognized option " << argv[0] << std::endl; + usage(opt.program); + return false; + } + } + if(argc > 0){ + std::cerr << "unrecognized option " << argv[0] << std::endl; + usage(opt.program); + return false; + } + return true; +} + +template +void reconstruction_test(PointIterator point_begin, PointIterator + point_end, TripleOutputIterator out, bool filter_input_points=false, + double perimeter=0) +{ + Options opt; + opt.abs_perimeter = perimeter; + std::cerr << "Compute Delaunay Tetrahedrization " << std::endl; + CGAL::Timer t1; + t1.start(); + + Triangulation_3 dt( boost::make_transform_iterator(point_begin, Auto_count()), + boost::make_transform_iterator(point_end, Auto_count() ) ); + t1.stop(); + std::cerr << " Inserted " << dt.number_of_vertices() << " points, " + << dt.number_of_cells() << " cells computed in " + << t1.time() << " sec." << std::endl; + + t1.reset(); + t1.start(); + Surface S(dt, opt); + t1.stop(); + std::cerr << "Reconstruction takes " << t1.time() << " sec.\n"; + std::cerr << " " << S.number_of_outliers() << " outliers.\n"; + std::cerr << " Reconstructed surface: " << S.number_of_facets() << + " facets, " << S.number_of_vertices() << " vertices.\n"; + std::cerr << " " << S.number_of_border_edges() << + " border edges.\n"; + std::cerr << " number of connected components <= " + << (std::max)(1, S.number_of_connected_components()-1) + << std::endl; + write_triple_indices(out, S); +} + + +//___________________________________________ +int main(int argc, char* argv[]) +{ + CGAL::Timer timer, total; + total.start(); + timer.start(); + //parse command line + Options opt; + std::cerr << "Option line for this execution is :" << std::endl; + if (!parse(argc, argv, opt)) + exit(0); + std::cerr << std::endl << std::endl; + + std::vector points; + + file_input(opt, points); + +#if 1 + std::cerr << "Time for reading " << timer.time() << " sec." << std::endl; + std::vector > triples; + reconstruction_test(points.begin(), points.end(), std::back_inserter(triples)); + + + std::cout << triples.size() << std::endl; + for(int i = 0; i < triples.size(); ++i){ + std::cout << "3 " << triples[i].first << " " << triples[i].second << " " << triples[i].third << " " << std::endl; + } + +#else + + std::cerr << "Compute Delaunay Tetrahedrization " << std::endl; + CGAL::Timer t1; + t1.start(); + + Triangulation_3 dt( boost::make_transform_iterator(points.begin(),Auto_count()), + boost::make_transform_iterator(points.end(), Auto_count() ) ); + t1.stop(); + std::cerr << " Inserted " << dt.number_of_vertices() << " points, " + << dt.number_of_cells() << " cells computed in " + << t1.time() << " sec." << std::endl; + + if (dt.dimension() < 3) { + std::cerr << "-- 2D sample of points ???" << std::endl; + exit(0); + } + + points.clear(); + + + Surface S(dt, opt); + + std::cerr << "Total time: " << timer.time() << " sec." << std::endl; + // write_to_file_vrml2(opt.foutname, S, opt.contour, opt.red, opt.green, opt.blue, opt.no_header); + write_to_file(opt.foutname, S, opt.contour, opt.out_format, opt.red, opt.green, opt.blue, opt.no_header); + + std::cerr << " " << S.number_of_outliers() + << " outliers." << std::endl; + std::cerr << " Reconstructed surface: " << S.number_of_facets() << + " facets, " << S.number_of_vertices() << " vertices." << std::endl; + std::cerr << " " << S.number_of_border_edges() << + " border edges." << std::endl; + std::cerr << " number of connected components <= " + << (std::max)(1, S.number_of_connected_components()-1) + << std::endl << std::endl; + +#endif + total.stop(); + std::cerr << "Total = " << total.time() << " sec." << std::endl; + return 0; +} + + + + diff --git a/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/fill.cpp b/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/fill.cpp new file mode 100644 index 00000000000..91b210b2f82 --- /dev/null +++ b/Advancing_front_surface_reconstruction/examples/Advancing_front_surface_reconstruction/fill.cpp @@ -0,0 +1,138 @@ +#include +#include +#include +#include +#include +#include +#include + + +typedef double NT; + +struct K : public CGAL::Filtered_kernel > {}; +typedef K::Point_3 Point; +typedef K::Vector_3 Vector; +typedef K::Segment_3 Segment; +typedef K::Triangle_3 Triangle; + +NT +weight(const Point& p, const Point& q, const Point& r){ + NT area = std::sqrt(Triangle(p,q,r).squared_area()); + NT l1 = std::sqrt((p-q) * (p-q)); + NT l2 = std::sqrt((p-r) * (p-r)); + NT l3 = std::sqrt((q-r) * (q-r)); + if(l1>l2) std::swap(l1,l2); + if(l2>l3) std::swap(l2,l3); + if(l1>l2) std::swap(l1,l2); + if(l2>l3) std::swap(l2,l3); + + // Taken from Piecewise-Linear Interpolation between Polygonal Slices + // from Gill Barequet and Micha Sharir + return 0.85 * area + 0.05 * (l1 +l2 + l3) + 0.1 * l3/l1; +} + + +void +insert(std::set >& triangles, int i, int j, int k){ + std::cout << i << ", " << j << ", " << k << std::endl; + if(i>j) std::swap(i,j); + if(j>k) std::swap(j,k); + if(i>j) std::swap(i,j); + if(j>k) std::swap(j,k); + std::cout << i << ", " << j << ", " << k << std::endl; + triangles.insert(CGAL::make_triple(i,j,k)); +} + + +void +collect(int i, int k, int n, const std::vector& O, std::set >& triangles){ + + std::cout << "collect(" << i << ", " << k << ")" << std::endl; + if((i+2) == k){ + insert(triangles, i, i+1, k); + }else { + int o = O[i*n+k]; + + if(o != (i+1)){ + collect(i, o, n, O, triangles); + } + insert(triangles, i, o, k); + if(o != (k-1)){ + collect(o, k, n, O, triangles); + } + } +} + + + +int +main(){ + + int n; + std::cin >> n; + + std::vector points(n); + CGAL::copy_n(std::istream_iterator(std::cin), n, points.begin()); + + + std::set > triangles; + + std::vector W(n*n); + std::vector O(n*n); + + for(int i = 0; i <= n-2; i++){ + W[i*n + i + 1] = 0; + } + for(int i = 0; i <= n-3; i++){ + W[i*n + i + 2] = weight(points[i], points[i+1], points[i+3]); + } + + for(int j = 3; j <= n-1; j++){ + for(int i=0; i <= n - j - 1; i++){ + int k = i + j; + double lmin = -1; + int lmin_index; + for(int m = i+1; m < k; m++){ + double d = W[i*n + m] + W[m*n + k] + weight(points[i], points[m], points[k]); + if( (lmin == -1) || (d < lmin )){ + lmin = d; + lmin_index = m; + } + } + W[i*n + k] = lmin; + O[i*n + k] = lmin_index; + + } + } + + collect(0, n-1, n, O, triangles); + + std::cout << "Shape {\n" + "appearance Appearance {\n" + "material Material { diffuseColor .9 .5 .1}}\n" + "geometry\n" + "IndexedFaceSet {\n" + "coord DEF def_coords Coordinate {\n" + "point [ \n" ; + + for (int i = 0; i < n; i++){ + std::cout << points[i].x() << " " << points[i].y() << " " << points[i].z() << ",\n "; + } + std::cout << "]\n" + "}\n" + "solid FALSE\n" + "coordIndex [ \n"; + + for(std::set >::iterator it = triangles.begin(); + it != triangles.end(); + it++){ + std::cout << it->first << ", " << it->second << ", " << it->third << ", -1,\n "; + } + std::cout << "]\n" + "}# IndexedFaceSet\n" + "}# Shape\n"; + + return 1; +} + + diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR/Surface_face_base_2.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/Surface_face_base_2.h new file mode 100644 index 00000000000..d59bf2048c6 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/Surface_face_base_2.h @@ -0,0 +1,88 @@ +#ifndef CGAL_AFSR_SURFACE_FACE_BASE_2_H +#define CGAL_AFSR_SURFACE_FACE_BASE_2_H + +// This face class stores a facet in the tetrahedrization +// When it gets reoriented by the TDS, it also changes the facet + +namespace CGAL { +namespace AFSR { + +template < typename GT, + typename F3, + typename Fb = CGAL::Triangulation_ds_face_base_2<> > +class Surface_face_base_2 + : public Fb +{ + typedef typename Fb::Triangulation_data_structure Tds; + +public: + typedef typename Tds::Face_handle Face_handle; + typedef typename Tds::Vertex_handle Vertex_handle; + + template < typename TDS2 > + struct Rebind_TDS { + typedef typename Fb::template Rebind_TDS::Other Fb2; + typedef Surface_face_base_2 Other; + }; + +private: + F3 _facet; + bool _is_on_surface; + +public: + Surface_face_base_2() + : Fb(), _is_on_surface(true) + {} + + Surface_face_base_2(Vertex_handle v0, + Vertex_handle v1, + Vertex_handle v2) + : Fb(v0, v1, v2), _is_on_surface(true) + {} + + Surface_face_base_2(Vertex_handle v0, + Vertex_handle v1, + Vertex_handle v2, + Face_handle n0, + Face_handle n1, + Face_handle n2) + : Fb(v0, v1, v2, n0, n1, n2), _is_on_surface(true) + {} + + void set_facet(const F3& facet) + { + _facet = facet; + } + + const F3& facet() const + { + return _facet; + } + + void set_is_on_surface(bool is_on_surface) + { + _is_on_surface = is_on_surface; + } + + bool is_on_surface() const + { + return _is_on_surface; + } + + + void reorient() + { + Fb::reorient(); + if( is_on_surface()){ + _facet = std::make_pair(_facet.first->neighbor(_facet.second), + _facet.first->neighbor(_facet.second)->index(_facet.first)); + } + } + +}; + + +} // namespace AFSR +} // namespace CGAL + +#endif // CGAL_AFSR_SURFACE_FACE_BASE_2_H diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR/Surface_vertex_base_2.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/Surface_vertex_base_2.h new file mode 100644 index 00000000000..d274a07cef7 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/Surface_vertex_base_2.h @@ -0,0 +1,76 @@ +// Copyright (c) 2005 GeometryFactory Sarl +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org); you may redistribute it under +// the terms of the Q Public License version 1.0. +// See the file LICENSE.QPL distributed with CGAL. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// +// Author(s) : Andreas Fabri + + +#ifndef CGAL_AFSR_SURFACE_VERTEX_BASE_2_H +#define CGAL_AFSR_SURFACE_VERTEX_BASE_2_H + +#include +#include + +namespace CGAL { +namespace AFSR { + +template < typename GT, + typename V3, + typename Vb = CGAL::Triangulation_ds_vertex_base_2<> > +class Surface_vertex_base_2 + : public Vb + +{ + typedef typename Vb::Triangulation_data_structure Tds; +public: + typedef GT Geom_traits; + typedef typename GT::Point_3 Point; + typedef Tds Triangulation_data_structure; + typedef typename Tds::Face_handle Face_handle; + typedef typename Tds::Vertex_handle Vertex_handle; + + template < typename TDS2 > + struct Rebind_TDS { + typedef typename Vb::template Rebind_TDS::Other Vb2; + typedef Surface_vertex_base_2 Other; + }; + +private: + V3 _vertex; +public: + Surface_vertex_base_2() : Vb() {} + Surface_vertex_base_2(Face_handle f) : Vb(f) {} + + void set_vertex(const V3& v) + { + _vertex = v; + } + + V3 vertex_3() const + { + return _vertex; + } + + const Point& point() const { return _vertex->point(); } + + +}; + + + + + +} // namespace AFSR +} // namespace CGAL + +#endif //CGAL::AFSR_SURFACE_VERTEX_BASE_2_H diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR/construct_surface_2.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/construct_surface_2.h new file mode 100644 index 00000000000..6ee1003ceb7 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/construct_surface_2.h @@ -0,0 +1,127 @@ +#ifndef CGAL_AFSR_CONSTRUCT_SURFACE_2 + +namespace CGAL { + +template +class Advancing_front_surface_reconstruction; + +namespace AFSR { + + +template +typename TDS::Vertex_handle +construct_surface(TDS& tds, const CGAL::Advancing_front_surface_reconstruction& surface) +{ + + typedef typename TDS::Vertex_handle Vertex_handle; + typedef std::pair Vh_pair; + typedef typename TDS::Face_handle Face_handle; + typedef typename TDS::Edge Edge; + + + Triangulation& T = surface.triangulation(); + // create an infinite-vertex and infinite faces with the + // boundary edges if any. + // return the infinite vertex if created + Vertex_handle vinf; + + std::vector vvh; + if(tds.number_of_vertices() != 0){ + tds.clear(); + } + int dim = 2; + tds.set_dimension(dim); + + CGAL::Unique_hash_map vertex_index_map(-1, T.number_of_vertices()); + + int i=0; + for (typename Triangulation::Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++){ + typename CGAL::Unique_hash_map::Data& d = vertex_index_map[v_it]; + if ((!v_it->is_exterior()) && d == -1){ + d = i; + Vertex_handle vh = tds.create_vertex(); + vvh.push_back(vh); + vh->set_vertex(typename Triangulation::Vertex_handle(v_it)); + i++; + } + } + std::map edge_map; + + for(typename Triangulation::Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + typename Triangulation::Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + n = c->neighbor(ci); + ni = n->index(c); + int i1, i2 ,i3; + + if (c->is_selected_facet(ci)) + { + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + + Face_handle fh = tds.create_face(vvh[vertex_index_map[c->vertex(i1)]], + vvh[vertex_index_map[c->vertex(i2)]], + vvh[vertex_index_map[c->vertex(i3)]]); + fh->set_facet(*f_it); + vvh[vertex_index_map[c->vertex(i1)]]->set_face(fh); + vvh[vertex_index_map[c->vertex(i2)]]->set_face(fh); + vvh[vertex_index_map[c->vertex(i3)]]->set_face(fh); + for (int ih = 0; ih < 3; ++ih) { + tds.set_adjacency(fh, ih, edge_map); + } + } + + if (n->is_selected_facet(ni)) + { + i1 = (ni+1) & 3; + i2 = (ni+2) & 3; + i3 = (ni+3) & 3; + Face_handle fh = tds.create_face(vvh[vertex_index_map[n->vertex(i1)]], + vvh[vertex_index_map[n->vertex(i2)]], + vvh[vertex_index_map[n->vertex(i3)]]); + fh->set_facet(std::make_pair(n, ni)); + vvh[vertex_index_map[n->vertex(i1)]]->set_face(fh); + vvh[vertex_index_map[n->vertex(i2)]]->set_face(fh); + vvh[vertex_index_map[n->vertex(i3)]]->set_face(fh); + for (int ih = 0; ih < 3; ++ih) { + tds.set_adjacency(fh, ih, edge_map); + } + } + + } + + if ( !edge_map.empty()) { + vinf = tds.create_vertex(); + vinf->set_vertex(T.infinite_vertex()); + std::map inf_edge_map; + while (!edge_map.empty()) { + Face_handle fh = edge_map.begin()->second.first; + int ih = edge_map.begin()->second.second; + Face_handle fn = tds.create_face( vinf, + fh->vertex(TDS::cw(ih)), + fh->vertex(TDS::ccw(ih))); + fn->set_facet(std::make_pair( typename Triangulation::Cell_handle(),0)); + fn->set_is_on_surface(false); + vinf->set_face(fn); + tds.set_adjacency(fn, 0, fh, ih); + tds.set_adjacency(fn, 1, inf_edge_map); + tds.set_adjacency(fn, 2, inf_edge_map); + edge_map.erase(edge_map.begin()); + } + CGAL_triangulation_assertion(inf_edge_map.empty()); + } + tds.reorient_faces(); + return vinf; + +} + + +} // namespace AFSR +} // namespace CGAL +#endif // CGAL_AFSR_CONSTRUCT_SURFACE_2 diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR/orient.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/orient.h new file mode 100644 index 00000000000..c74232afe88 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR/orient.h @@ -0,0 +1,120 @@ +#ifndef CGAL_AFSR_ORIENT + +namespace CGAL { +namespace AFSR { + + +template +typename TDS::Vertex_handle +orient(TDS& tds, const Advancing_front_surface_reconstruction& surface) +{ + + typedef typename TDS::Vertex_handle Vertex_handle; + typedef std::pair Vh_pair; + typedef typename TDS::Face_handle Face_handle; + typedef typename TDS::Edge Edge; + + Triangulation& T = surface.triangulation(); + // create an infinite-vertex and infinite faces with the + // boundary edges if any. + // return the infinite vertex if created + Vertex_handle vinf; + + std::vector vvh; + if(tds.number_of_vertices() != 0) tds.clear(); + int dim = 2; + tds.set_dimension(dim); + + CGAL::Unique_hash_map vertex_index_map(-1, T.number_of_vertices()); + + int i=0; + for (typename Triangulation::Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++){ + typename CGAL::Unique_hash_map::Data& d = vertex_index_map[v_it]; + if ((!v_it->is_exterior()) && d == -1){ + d = i; + Vertex_handle vh = tds.create_vertex(); + vvh.push_back(vh); + vh->set_point(v_it->point()); + i++; + } + } + std::map edge_map; + + + for(typename Triangulation::Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + typename Triangulation::Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + n = c->neighbor(ci); + ni = n->index(c); + int i1, i2 ,i3; + + if (c->is_selected_facet(ci)) + { + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + + Face_handle fh = tds.create_face(vvh[vertex_index_map[c->vertex(i1)]], + vvh[vertex_index_map[c->vertex(i2)]], + vvh[vertex_index_map[c->vertex(i3)]]); + vvh[vertex_index_map[c->vertex(i1)]]->set_face(fh); + vvh[vertex_index_map[c->vertex(i2)]]->set_face(fh); + vvh[vertex_index_map[c->vertex(i3)]]->set_face(fh); + for (int ih = 0; ih < 3; ++ih) { + tds.set_adjacency(fh, ih, edge_map); + } + } + + if (n->is_selected_facet(ni)) + { + i1 = (ni+1) & 3; + i2 = (ni+2) & 3; + i3 = (ni+3) & 3; + Face_handle fh = tds.create_face(vvh[vertex_index_map[n->vertex(i1)]], + vvh[vertex_index_map[n->vertex(i2)]], + vvh[vertex_index_map[n->vertex(i3)]]); + vvh[vertex_index_map[n->vertex(i1)]]->set_face(fh); + vvh[vertex_index_map[n->vertex(i2)]]->set_face(fh); + vvh[vertex_index_map[n->vertex(i3)]]->set_face(fh); + for (int ih = 0; ih < 3; ++ih) { + tds.set_adjacency(fh, ih, edge_map); + } + } + + } + + if ( !edge_map.empty()) { + vinf = tds.create_vertex(); + std::map inf_edge_map; + while (!edge_map.empty()) { + Face_handle fh = edge_map.begin()->second.first; + int ih = edge_map.begin()->second.second; + Face_handle fn = tds.create_face( vinf, + fh->vertex(TDS::cw(ih)), + fh->vertex(TDS::ccw(ih))); + vinf->set_face(fn); + tds.set_adjacency(fn, 0, fh, ih); + tds.set_adjacency(fn, 1, inf_edge_map); + tds.set_adjacency(fn, 2, inf_edge_map); + edge_map.erase(edge_map.begin()); + } + CGAL_triangulation_assertion(inf_edge_map.empty()); + } + + + tds.reorient_faces(); + return vinf; + +} + + + +} // namespace AFSR +} // namespace CGAL + +#endif //CGAL_AFSR_ORIENT diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR_cell_base_3.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_cell_base_3.h new file mode 100644 index 00000000000..0c279cda67c --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_cell_base_3.h @@ -0,0 +1,224 @@ +#ifndef CGAL_AFSR_CELL_BASE_3_H +#define CGAL_AFSR_CELL_BASE_3_H + +namespace CGAL { + +template < class CellBase > +class AFSR_cell_base_3 : public CellBase +{ + +public: + template < typename TDS2 > + struct Rebind_TDS { + typedef typename CellBase::template Rebind_TDS::Other Cb2; + typedef AFSR_cell_base_3 Other; + }; + + typedef typename CellBase::Vertex_handle Vertex_handle; + typedef typename CellBase::Cell_handle Cell_handle; + +private: + +#ifdef AFSR_FACET_NUMBER + int _facet_number[4]; +#endif + typedef double coord_type; +#ifdef AFSR_LAZY + typedef typename CGAL::Simple_cartesian::Point_3 D_Point; +#endif + //-------------------- DATA MEMBERS --------------------------------- + + coord_type* _smallest_radius_facet_tab; + unsigned char selected_facet; +#ifdef AFSR_LAZY + D_Point* _circumcenter; + coord_type* _squared_radius; +#endif + + //-------------------- CONSTRUCTORS ---------------------------------- + +public: + + AFSR_cell_base_3() + : CellBase(), + _smallest_radius_facet_tab(NULL), selected_facet(0) +#ifdef AFSR_LAZY + , _circumcenter(NULL), _squared_radius(NULL) +#endif + { +#ifdef AFSR_FACET_NUMBER + for(int i = 0; i < 4; i++){ + _facet_number[i] = -1; + } +#endif + } + + AFSR_cell_base_3(Vertex_handle v0, Vertex_handle v1, Vertex_handle v2, Vertex_handle v3) + : CellBase( v0, v1, v2, v3), + _smallest_radius_facet_tab(NULL), selected_facet(0) +#ifdef AFSR_LAZY + , _circumcenter(NULL), _squared_radius(NULL) +#endif + { +#ifdef FACET_NUMBER + for(int i = 0; i < 4; i++){ + _facet_number[i] = -1; + } +#endif + } + + AFSR_cell_base_3(Vertex_handle v0, Vertex_handle v1, Vertex_handle v2, Vertex_handle v3, + Cell_handle n0, Cell_handle n1, Cell_handle n2, Cell_handle n3) + : CellBase(v0, v1, v2, v3, + n0, n1, n2, n3), + _smallest_radius_facet_tab(NULL), selected_facet(0) +#ifdef AFSR_LAZY + , _circumcenter(NULL), _squared_radius(NULL) +#endif + { +#ifdef AFSR_FACET_NUMBER + for(int i = 0; i < 4; i++){ + _facet_number[i] = -1; + } +#endif + } + + //-------------------- DESTRUCTOR ----------------------------------- + + inline ~AFSR_cell_base_3() + { + if (_smallest_radius_facet_tab != NULL) + delete[] _smallest_radius_facet_tab; +#ifdef AFSR_LAZY + if (_circumcenter != NULL) + delete _circumcenter; + if (_squared_radius != NULL) + delete _squared_radius; +#endif + } + + //-------------------- MEMBER FUNCTIONS ---------------------------- + +public: + + inline void clear() + { + if (_smallest_radius_facet_tab != NULL) + delete[] _smallest_radius_facet_tab; + _smallest_radius_facet_tab = NULL; + selected_facet = 0; +#ifdef AFSR_LAZY + if (_circumcenter != NULL) + delete _circumcenter; + _circumcenter = NULL; + if (_squared_radius != NULL) + delete _squared_radius; + _squared_radius = NULL; +#endif + } + + //------------------------------------------------------------------- + + inline coord_type get_smallest_radius(const int& i) + { + if (_smallest_radius_facet_tab == NULL) + return -1; + return _smallest_radius_facet_tab[i]; + } + + inline void set_smallest_radius(const int& i, const coord_type& c) + { + if (_smallest_radius_facet_tab == NULL) + { + _smallest_radius_facet_tab = new coord_type[4]; + for(int i = 0; i < 4; i++) + _smallest_radius_facet_tab[i] = -1; + } + _smallest_radius_facet_tab[i] = c; + } + + // pour un controle de l'allocation memoire... utile??? + inline bool alloc_smallest_radius_tab(coord_type* ptr) + { + if (_smallest_radius_facet_tab==NULL) + { + _smallest_radius_facet_tab = ptr; + return true; + } + return false; + } + + + //------------------------------------------------------------------- + +#ifdef FACET_NUMBER + void set_facet_number(int i, int n){} + { + _facet_number[i] = n; + } + + int facet_number(int i) + { + return _facet_number[i]; + } +#else + void set_facet_number(int, int){} + int facet_number(int){return 0;} +#endif + + + //------------------------------------------------------------------- + + inline void select_facet(const int& i) + { + selected_facet |= (1 << i); + } + + inline void unselect_facet(const int& i) + { + selected_facet &= (15 - (1 << i)); + } + + inline bool is_selected_facet(const int& i) + { + return (selected_facet & (1 << i)) != 0; + } + + inline bool has_facet_on_surface(const int& i) + { + return (selected_facet & (1 << i)) != 0; + } + +#ifdef AFSR_LAZY + + //------------------------------------------------------------------- + + inline D_Point* get_lazy_circumcenter() + { + return _circumcenter; + } + + inline void set_lazy_circumcenter(const D_Point& p) + { + _circumcenter = new D_Point(p); + } + + //------------------------------------------------------------------- + + inline coord_type* get_lazy_squared_radius() + { + return _squared_radius; + } + + inline void set_lazy_squared_radius(const coord_type& sr) + { + _squared_radius = new coord_type(sr); + } + +#endif //AFSR_LAZY + +}; + +} // namespace CGAL + +#endif // CGAL_AFSR_CELL_BASE_3_H diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR_options.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_options.h new file mode 100644 index 00000000000..3fda69a12c2 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_options.h @@ -0,0 +1,46 @@ +#ifndef CGAL_AFSR_OPTIONS_H +#define CGAL_AFSR_OPTIONS_H + + +namespace CGAL { + +class AFSR_options { +public: + AFSR_options() + : file_input(true), file_output(false), + Delaunay(false), contour(false), binary(false), xyz(false), + Section_file(false), max_connected_comp(-1), + delta(.86), K_init(1.1), K_step(.1), K(5), out_format(0), + NB_BORDER_MAX(15), red(0), green(0), blue(0), no_header(false), area(0), perimeter(0), + abs_area(0), abs_perimeter(0) + { + std::strcpy(finname,"finput"); + std::strcpy(foutname,"foutput"); + } + + char program[100]; + char finname[100]; + char foutname[100]; + bool file_input; + bool file_output; + bool Delaunay; + bool contour; + bool binary; + bool xyz; + bool Section_file; + int max_connected_comp; + double delta; + double K_init; + double K_step; + double K; + int out_format; + int NB_BORDER_MAX; + double red, green, blue; + bool no_header; + double area, perimeter, abs_area, abs_perimeter; +}; + + +} // namespace CGAL + +#endif // CGAL_AFSR_OPTIONS_H diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR_vertex_base_3.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_vertex_base_3.h new file mode 100644 index 00000000000..1820698a372 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_vertex_base_3.h @@ -0,0 +1,473 @@ +#ifndef CGAL_AFSR_VERTEX_BASE_3_H +#define CGAL_AFSR_VERTEX_BASE_3_H + +#include + +#include +#include +#include + +#include + +namespace CGAL { + +template > +class AFSR_vertex_base_3 : public VertexBase +{ +public: + + template < typename TDS2 > + struct Rebind_TDS { + typedef typename VertexBase::template Rebind_TDS::Other Vb2; + typedef AFSR_vertex_base_3 Other; + }; + + + typedef VertexBase Base; + typedef typename Base::Vertex_handle Vertex_handle; + typedef typename Base::Cell_handle Cell_handle; + typedef typename VertexBase::Point Point; + typedef double coord_type; + + typedef Triple< Cell_handle, int, int > Edge; + typedef std::pair< Edge, int > Edge_incident_facet; + typedef std::pair< Edge_incident_facet, Edge_incident_facet > IO_edge_type; + + typedef coord_type criteria; + + typedef std::pair< criteria, IO_edge_type > Radius_edge_type; + typedef std::pair< Radius_edge_type, int > Border_elt; + typedef std::pair< Vertex_handle, Border_elt > Next_border_elt; + +private: + + //par convention je remplis d'abord first et si necessaire second... + typedef std::pair< Next_border_elt*, Next_border_elt*> Intern_successors_type; + +public: + + typedef std::pair< criteria, IO_edge_type* > Radius_ptr_type; + typedef std::pair< Vertex_handle, Vertex_handle > Edge_like; + typedef std::pair< criteria, Edge_like > Incidence_request_elt; + typedef std::list< Incidence_request_elt > Incidence_request_type; + typedef typename Incidence_request_type::iterator Incidence_request_iterator; + + // typedef std::set< void* > Interior_edge_set_type; + + //-------------------- DATA MEMBERS --------------------------------- + +private: + + int _mark; + int _post_mark; + Intern_successors_type* _incident_border; + + // Instead of having a set per vertex, there is a global list. + static std::list interior_edges; + // and two iterators per vertex in this list + // Note that ie_last is not past the end + // ie_first == ie_last == interior_edge.end() iff the set is empty + typename std::list::iterator ie_first, ie_last; + + + // We do the same for the incidence requests + static std::list< Incidence_request_elt > incidence_requests; + typename std::list< Incidence_request_elt >::iterator ir_first, ir_last; + //-------------------- CONSTRUCTORS --------------------------------- + +public: + + AFSR_vertex_base_3() + : VertexBase(), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new Intern_successors_type(new Next_border_elt(), + new Next_border_elt()); + _incident_border->first->first = NULL; + _incident_border->second->first = NULL; + } + + AFSR_vertex_base_3(const Point & p) + : VertexBase(p), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new Intern_successors_type(new Next_border_elt(), + new Next_border_elt()); + _incident_border->first->first = NULL; + _incident_border->second->first = NULL; + } + + AFSR_vertex_base_3(const Point & p, Cell_handle f) + : VertexBase(p, f), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new Intern_successors_type(new Next_border_elt(), + new Next_border_elt()); + _incident_border->first->first = NULL; + _incident_border->second->first = NULL; + } + + AFSR_vertex_base_3(Cell_handle f) + : VertexBase(f), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new Intern_successors_type(new Next_border_elt(), + new Next_border_elt()); + _incident_border->first->first = NULL; + _incident_border->second->first = NULL; + } + + AFSR_vertex_base_3(const AFSR_vertex_base_3& other) + : VertexBase(), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new Intern_successors_type(new Next_border_elt(), + new Next_border_elt()); + _incident_border->first->first = NULL; + _incident_border->second->first = NULL; + } + //-------------------- DESTRUCTOR ----------------------------------- + + ~AFSR_vertex_base_3() + { + if (_incident_border != NULL) + { + delete _incident_border->first; + delete _incident_border->second; + delete _incident_border; + } + if(ir_first != incidence_requests.end()){ + assert(ir_last != incidence_requests.end()); + typename std::list< Incidence_request_elt >::iterator b(ir_first), e(ir_last); + e++; + incidence_requests.erase(b, e); + } + + if(ie_first != interior_edges.end()){ + assert(ie_last != interior_edges.end()); + typename std::list::iterator b(ie_first), e(ie_last); + e++; + interior_edges.erase(b, e); + } + } + + //-------------------- MEMBER FUNCTIONS ----------------------------- + + inline void re_init() + { + if (_incident_border != NULL) + { + delete _incident_border->first; + delete _incident_border->second; + delete _incident_border; + } + + if(ir_first != incidence_requests.end()){ + assert(ir_last != incidence_requests.end()); + typename std::list< Incidence_request_elt >::iterator b(ir_first), e(ir_last); + e++; + incidence_requests.erase(b, e); + ir_first = incidence_requests.end(); + ir_last = incidence_requests.end(); + } + + _incident_border = new Intern_successors_type(new Next_border_elt(), + new Next_border_elt()); + _incident_border->first->first = NULL; + _incident_border->second->first = NULL; + _mark = -1; + _post_mark = -1; + } + + //------------------------------------------------------------------- + + inline bool is_on_border(const int& i) const + { + if (_incident_border == NULL) return false; //vh is interior + if (_incident_border->first->first != NULL) + { + if (_incident_border->second->first != NULL) + return ((_incident_border->first->second.second == i)|| + (_incident_border->second->second.second == i)); + return (_incident_border->first->second.second == i); + } + return false; //vh is still exterior + } + + inline Next_border_elt* get_next_on_border(const int& i) const + { + if (_incident_border == NULL) return NULL; //vh is interior + if (_incident_border->first->first != NULL) + if (_incident_border->first->second.second == i) + return _incident_border->first; + if (_incident_border->second->first != NULL) + if (_incident_border->second->second.second == i) + return _incident_border->second; + return NULL; + } + + + inline void remove_border_edge(Vertex_handle v) + { + if (_incident_border != NULL) + { + if (_incident_border->second->first == v) + { + _incident_border->second->first = NULL; + set_interior_edge(v); + return; + } + if (_incident_border->first->first == v) + { + if (_incident_border->second->first != NULL) + { + Next_border_elt* tmp = _incident_border->first; + _incident_border->first = _incident_border->second; + _incident_border->second = tmp; + _incident_border->second->first = NULL; + set_interior_edge(v); + return; + } + else + { + _incident_border->first->first = NULL; + set_interior_edge(v); + return; + } + } + } + } + + inline bool is_border_edge(Vertex_handle v) const + { + if (_incident_border == NULL) return false; + return ((_incident_border->first->first == v)|| + (_incident_border->second->first == v)); + } + + inline Next_border_elt* get_border_elt(Vertex_handle v) const + { + if (_incident_border == NULL) return NULL; + if (_incident_border->first->first == v) return _incident_border->first; + if (_incident_border->second->first == v) return _incident_border->second; + return NULL; + } + + inline Next_border_elt* first_incident() const + { + if (_incident_border == NULL) return NULL; + return _incident_border->first; + } + + inline Next_border_elt* second_incident() const + { + if (_incident_border == NULL) return NULL; + return _incident_border->second; + } + + + inline void set_next_border_elt(const Next_border_elt& elt) + { + if (_incident_border->first->first == NULL) + *_incident_border->first = elt; + else + { + if (_incident_border->second->first != NULL) + std::cerr << "+++probleme de MAJ du bord " << std::endl; + *_incident_border->second = elt; + } + } + + //------------------------------------------------------------------- + // pour gerer certaines aretes interieures: a savoir celle encore connectee au + // bord (en fait seule, les aretes interieures reliant 2 bords nous + // interressent...) + + inline bool is_interior_edge(Vertex_handle v) const + { + + bool r1; + if(ie_first == interior_edges.end()){ + r1 = false; + }else { + typename std::list::iterator b(ie_first), e(ie_last); + e++; + typename std::list::iterator r = std::find(b, e, v); + r1 = ( r != e); + } + + return r1; + } + + inline void set_interior_edge(Vertex_handle v) + { + if(ie_last == interior_edges.end()){ // empty set + assert(ie_first == ie_last); + ie_last = interior_edges.insert(ie_last, v); + ie_first = ie_last; + } else { + typename std::list::iterator e(ie_last); + e++; +#ifdef DEBUG + typename std::list::iterator r = std::find(ie_first, e, v); +#endif + assert(r == e); + ie_last = interior_edges.insert(e, v); + } + } + + inline void remove_interior_edge(Vertex_handle v) + { + if(ie_first == interior_edges.end()){ + assert(ie_last == ie_first); + } else if(ie_first == ie_last){ // there is only one element + if(*ie_first == v){ + interior_edges.erase(ie_first); + ie_last = interior_edges.end(); + ie_first = ie_last; + } + } else { + typename std::list::iterator b(ie_first), e(ie_last); + e++; + typename std::list::iterator r = std::find(b, e, v); + if(r != e){ + if(r == ie_first){ + ie_first++; + } + if(r == ie_last){ + ie_last--; + } + interior_edges.erase(r); + } + } + } + + //------------------------------------------------------------------- + + inline void set_incidence_request(const Incidence_request_elt& ir) + { + if(ir_last == incidence_requests.end()){ + assert(ir_first == ir_last); + ir_last = incidence_requests.insert(ir_last, ir); + ir_first = ir_last; + } else { + typename std::list::iterator e(ir_last); + e++; + ir_last = incidence_requests.insert(e, ir); + } + } + + inline bool is_incidence_requested() const + { + if(ir_last == incidence_requests.end()){ + assert(ir_first == incidence_requests.end()); + } + return (ir_last != incidence_requests.end()); + } + + inline Incidence_request_iterator incidence_request_begin() + { + return ir_first; + } + + inline Incidence_request_iterator get_incidence_request_end() + { + if(ir_last != incidence_requests.end()){ + assert(ir_first != incidence_requests.end()); + Incidence_request_iterator it(ir_last); + it++; + return it; + } + return ir_last; + } + + inline void erase_incidence_request() + { + if(ir_last != incidence_requests.end()){ + assert(ir_first != incidence_requests.end()); + ir_last++; + incidence_requests.erase(ir_first, ir_last); + ir_first = incidence_requests.end(); + ir_last = incidence_requests.end(); + } + } + + + //------------------------------------------------------------------- + + inline bool is_on_border() const + { + return (_mark > 0); + } + + inline bool not_interior() const + { + return (_mark != 0); + } + + inline int number_of_incident_border() const + { + return _mark; + } + + inline bool is_exterior() const + { + return (_mark < 0); + } + + //------------------------------------------------------------------- + + inline void inc_mark() + { + if (_mark==-1) + _mark=1; + else + _mark++; + } + + inline void dec_mark() + { + _mark--; + if(_mark == 0) + { + delete _incident_border->first; + delete _incident_border->second; + delete _incident_border; + _incident_border = NULL; + erase_incidence_request(); + } + } + + //------------------------------------------------------------------- + + inline void set_post_mark(const int& i) + { + _post_mark = i; + } + + inline bool is_post_marked(const int& i) + { + return (_post_mark == i); + } +}; + +template +std::list::Vertex_handle> AFSR_vertex_base_3::interior_edges; + +template +std::list::Incidence_request_elt> AFSR_vertex_base_3::incidence_requests; + +} // namespace CGAL + +#endif //CGAL_AFSR_VERTEX_BASE_3_H + diff --git a/Advancing_front_surface_reconstruction/include/CGAL/AFSR_vertex_base_with_id_3.h b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_vertex_base_with_id_3.h new file mode 100644 index 00000000000..f5f279f7c35 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/AFSR_vertex_base_with_id_3.h @@ -0,0 +1,557 @@ +#ifndef CGAL_AFSR_VERTEX_BASE_WITH_ID_3_H +#define CGAL_AFSR_VERTEX_BASE_WITH_ID_3_H + +#include + +#include + +#include +#include +#include + +#include + + +namespace CGAL { + +template > +class AFSR_vertex_base_with_id_3 : public VertexBase +{ +public: + + template < typename TDS2 > + struct Rebind_TDS { + typedef typename VertexBase::template Rebind_TDS::Other Vb2; + typedef AFSR_vertex_base_with_id_3 Other; + }; + + + + + typedef VertexBase Base; + typedef typename Base::Vertex_handle Vertex_handle; + typedef typename Base::Cell_handle Cell_handle; + typedef typename VertexBase::Point Point; + typedef double coord_type; + + typedef Triple< Cell_handle, int, int > Edge; + typedef std::pair< Edge, int > Edge_incident_facet; + typedef std::pair< Edge_incident_facet, Edge_incident_facet > IO_edge_type; + + typedef coord_type criteria; + + typedef std::pair< criteria, IO_edge_type > Radius_edge_type; + typedef std::pair< Radius_edge_type, int > Border_elt; + typedef std::pair< Vertex_handle, Border_elt > Next_border_elt; + +private: + + //par convention je remplis d'abord first et si necessaire second... + typedef std::pair< Next_border_elt*, Next_border_elt*> Intern_successors_type; + + + +public: + + typedef std::pair< criteria, IO_edge_type* > Radius_ptr_type; + typedef std::pair< Vertex_handle, Vertex_handle > Edge_like; + typedef std::pair< criteria, Edge_like > Incidence_request_elt; + typedef std::list< Incidence_request_elt > Incidence_request_type; + typedef typename Incidence_request_type::iterator Incidence_request_iterator; + + // typedef std::set< void* > Interior_edge_set_type; + + //-------------------- DATA MEMBERS --------------------------------- + + typedef int Info; // so that we are a model of TriangulationVertexBaseWithInfo_3 + +private: + int _id; + int _mark; + int _post_mark; + Intern_successors_type* _incident_border; + + // Instead of having a set per vertex, there is a global list. + static std::list interior_edges; + // and two iterators per vertex in this list + // Note that ie_last is not past the end + // ie_first == ie_last == interior_edge.end() iff the set is empty + typename std::list::iterator ie_first, ie_last; + + + // We do the same for the incidence requests + static std::list< Incidence_request_elt > incidence_requests; + typename std::list< Incidence_request_elt >::iterator ir_first, ir_last; + + static std::list bbb; + static boost::object_pool nbe_pool; + static boost::object_pool ist_pool; + + + + + //-------------------- CONSTRUCTORS --------------------------------- + + +public: + + Intern_successors_type* new_border() + { + Intern_successors_type* ret; + /* + std::allocator nbea; + std::allocator ista; + Next_border_elt nbe; + Next_border_elt* p1 = nbea.allocate(1); + Next_border_elt* p2 = nbea.allocate(1); + nbea.construct(p1, nbe); + nbea.construct(p1, nbe); + + ret = ista.allocate(1); + ista.construct(ret, std::make_pair(p1, p2)); + */ + + Next_border_elt* p1 = nbe_pool.malloc(); + Next_border_elt* p2 = nbe_pool.malloc(); + ret = ist_pool.malloc(); + + Intern_successors_type ist(p1,p2); + *ret = ist; + + + /* + ret = new Intern_successors_type(new Next_border_elt(), + new Next_border_elt()); + */ + ret->first->first = NULL; + ret->second->first = NULL; + return ret; + } + + + void delete_border() + { + /* + std::allocator nbea; + std::allocator ista; + nbea.destroy(_incident_border->first); + nbea.destroy(_incident_border->second); + nbea.deallocate(_incident_border->first, 1); + nbea.deallocate(_incident_border->second, 1); + ista.destroy(_incident_border); + ista.deallocate(_incident_border, 1); + */ + /* + delete _incident_border->first; + delete _incident_border->second; + */ + //nbe_pool.free(_incident_border->first); + //nbe_pool.free(_incident_border->second); + //delete _incident_border; + + _incident_border = NULL; + + + } + + + AFSR_vertex_base_with_id_3() + : VertexBase(), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new_border(); + } + + AFSR_vertex_base_with_id_3(const Point & p) + : VertexBase(p), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new_border(); + } + + AFSR_vertex_base_with_id_3(const Point & p, Cell_handle f) + : VertexBase(p, f), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new_border(); + } + + AFSR_vertex_base_with_id_3(Cell_handle f) + : VertexBase(f), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new_border(); + } + + AFSR_vertex_base_with_id_3(const AFSR_vertex_base_with_id_3& other) + : VertexBase(), _mark(-1), + _post_mark(-1), + ie_first(interior_edges.end()), ie_last(interior_edges.end()), + ir_first(incidence_requests.end()), ir_last(incidence_requests.end()) + { + _incident_border = new_border(); + } + //-------------------- DESTRUCTOR ----------------------------------- + + ~AFSR_vertex_base_with_id_3() + { + if (_incident_border != NULL) + { + delete_border(); + } + if(ir_first != incidence_requests.end()){ + assert(ir_last != incidence_requests.end()); + typename std::list< Incidence_request_elt >::iterator b(ir_first), e(ir_last); + e++; + incidence_requests.erase(b, e); + } + + if(ie_first != interior_edges.end()){ + assert(ie_last != interior_edges.end()); + typename std::list::iterator b(ie_first), e(ie_last); + e++; + interior_edges.erase(b, e); + } + } + + //-------------------- MEMBER FUNCTIONS ----------------------------- + + int& id() + { + return _id; + } + + const int& id() const + { + return _id; + } + + int& info() + { + return _id; + } + + const int& info() const + { + return _id; + } + inline void re_init() + { + if (_incident_border != NULL) + { + delete_border(); + //delete _incident_border->first; + //delete _incident_border->second; + //delete _incident_border; + } + + if(ir_first != incidence_requests.end()){ + assert(ir_last != incidence_requests.end()); + typename std::list< Incidence_request_elt >::iterator b(ir_first), e(ir_last); + e++; + incidence_requests.erase(b, e); + ir_first = incidence_requests.end(); + ir_last = incidence_requests.end(); + } + + _incident_border = new_border(); + _mark = -1; + _post_mark = -1; + } + + //------------------------------------------------------------------- + + inline bool is_on_border(const int& i) const + { + if (_incident_border == NULL) return false; //vh is interior + if (_incident_border->first->first != NULL) + { + if (_incident_border->second->first != NULL) + return ((_incident_border->first->second.second == i)|| + (_incident_border->second->second.second == i)); + return (_incident_border->first->second.second == i); + } + return false; //vh is still exterior + } + + inline Next_border_elt* get_next_on_border(const int& i) const + { + if (_incident_border == NULL) return NULL; //vh is interior + if (_incident_border->first->first != NULL) + if (_incident_border->first->second.second == i) + return _incident_border->first; + if (_incident_border->second->first != NULL) + if (_incident_border->second->second.second == i) + return _incident_border->second; + return NULL; + } + + + inline void remove_border_edge(Vertex_handle v) + { + if (_incident_border != NULL) + { + if (_incident_border->second->first == v) + { + _incident_border->second->first = NULL; + set_interior_edge(v); + return; + } + if (_incident_border->first->first == v) + { + if (_incident_border->second->first != NULL) + { + Next_border_elt* tmp = _incident_border->first; + _incident_border->first = _incident_border->second; + _incident_border->second = tmp; + _incident_border->second->first = NULL; + set_interior_edge(v); + return; + } + else + { + _incident_border->first->first = NULL; + set_interior_edge(v); + return; + } + } + } + } + + inline bool is_border_edge(Vertex_handle v) const + { + if (_incident_border == NULL) return false; + return ((_incident_border->first->first == v)|| + (_incident_border->second->first == v)); + } + + inline Next_border_elt* get_border_elt(Vertex_handle v) const + { + if (_incident_border == NULL) return NULL; + if (_incident_border->first->first == v) return _incident_border->first; + if (_incident_border->second->first == v) return _incident_border->second; + return NULL; + } + + inline Next_border_elt* first_incident() const + { + if (_incident_border == NULL) return NULL; + return _incident_border->first; + } + + inline Next_border_elt* second_incident() const + { + if (_incident_border == NULL) return NULL; + return _incident_border->second; + } + + + inline void set_next_border_elt(const Next_border_elt& elt) + { + if (_incident_border->first->first == NULL) + *_incident_border->first = elt; + else + { + if (_incident_border->second->first != NULL) + std::cerr << "+++probleme de MAJ du bord " << std::endl; + *_incident_border->second = elt; + } + } + + //------------------------------------------------------------------- + // pour gerer certaines aretes interieures: a savoir celle encore connectee au + // bord (en fait seule, les aretes interieures reliant 2 bords nous + // interressent...) + + inline bool is_interior_edge(Vertex_handle v) const + { + + bool r1; + if(ie_first == interior_edges.end()){ + r1 = false; + }else { + typename std::list::iterator b(ie_first), e(ie_last); + e++; + typename std::list::iterator r = std::find(b, e, v); + r1 = ( r != e); + } + + return r1; + } + + inline void set_interior_edge(Vertex_handle v) + { + if(ie_last == interior_edges.end()){ // empty set + assert(ie_first == ie_last); + ie_last = interior_edges.insert(ie_last, v); + ie_first = ie_last; + } else { + typename std::list::iterator e(ie_last); + e++; +#ifdef DEBUG + typename std::list::iterator r = std::find(ie_first, e, v); + assert(r == e); +#endif + ie_last = interior_edges.insert(e, v); + } + } + + inline void remove_interior_edge(Vertex_handle v) + { + if(ie_first == interior_edges.end()){ + assert(ie_last == ie_first); + } else if(ie_first == ie_last){ // there is only one element + if(*ie_first == v){ + interior_edges.erase(ie_first); + ie_last = interior_edges.end(); + ie_first = ie_last; + } + } else { + typename std::list::iterator b(ie_first), e(ie_last); + e++; + typename std::list::iterator r = std::find(b, e, v); + if(r != e){ + if(r == ie_first){ + ie_first++; + } + if(r == ie_last){ + ie_last--; + } + interior_edges.erase(r); + } + } + } + + //------------------------------------------------------------------- + + inline void set_incidence_request(const Incidence_request_elt& ir) + { + if(ir_last == incidence_requests.end()){ + assert(ir_first == ir_last); + ir_last = incidence_requests.insert(ir_last, ir); + ir_first = ir_last; + } else { + typename std::list::iterator e(ir_last); + e++; + ir_last = incidence_requests.insert(e, ir); + } + } + + inline bool is_incidence_requested() const + { + if(ir_last == incidence_requests.end()){ + assert(ir_first == incidence_requests.end()); + } + return (ir_last != incidence_requests.end()); + } + + inline Incidence_request_iterator incidence_request_begin() + { + return ir_first; + } + + inline Incidence_request_iterator get_incidence_request_end() + { + if(ir_last != incidence_requests.end()){ + assert(ir_first != incidence_requests.end()); + Incidence_request_iterator it(ir_last); + it++; + return it; + } + return ir_last; + } + + inline void erase_incidence_request() + { + if(ir_last != incidence_requests.end()){ + assert(ir_first != incidence_requests.end()); + ir_last++; + incidence_requests.erase(ir_first, ir_last); + ir_first = incidence_requests.end(); + ir_last = incidence_requests.end(); + } + } + + + //------------------------------------------------------------------- + + inline bool is_on_border() const + { + return (_mark > 0); + } + + inline bool not_interior() const + { + return (_mark != 0); + } + + inline int number_of_incident_border() const + { + return _mark; + } + + inline bool is_exterior() const + { + return (_mark < 0); + } + + //------------------------------------------------------------------- + + inline void inc_mark() + { + if (_mark==-1) + _mark=1; + else + _mark++; + } + + inline void dec_mark() + { + _mark--; + if(_mark == 0) + { + delete_border(); + erase_incidence_request(); + } + } + + //------------------------------------------------------------------- + + inline void set_post_mark(const int& i) + { + _post_mark = i; + } + + inline bool is_post_marked(const int& i) + { + return (_post_mark == i); + } +}; + + +template +boost::object_pool::Next_border_elt> AFSR_vertex_base_with_id_3::nbe_pool; + +template +boost::object_pool::Intern_successors_type> AFSR_vertex_base_with_id_3::ist_pool; + + + +template +std::list::Vertex_handle> AFSR_vertex_base_with_id_3::interior_edges; + +template +std::list::Incidence_request_elt> AFSR_vertex_base_with_id_3::incidence_requests; + +} // namespace CGAL + +#endif // CGALAFSR_VERTEX_BASE_3_H + diff --git a/Advancing_front_surface_reconstruction/include/CGAL/Advancing_front_surface_reconstruction.h b/Advancing_front_surface_reconstruction/include/CGAL/Advancing_front_surface_reconstruction.h new file mode 100644 index 00000000000..af2a490c356 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/Advancing_front_surface_reconstruction.h @@ -0,0 +1,1955 @@ +#ifndef CGAL_ADVANCING_FRONT_SURFACE_RECONSTRUCTION_H +#define CGAL_ADVANCING_FRONT_SURFACE_RECONSTRUCTION_H + +// In order to activate lazy evaluation: +// #define LAZY + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace CGAL { + +// This iterator allows to visit all contours. It has the particularity +// that it visits the entry point of the contour twice. This allows to +// detect that the traversal of the border is finished. One more increment +// brings us to the next vertex. + + + +template < class Surface> +class Advancing_front_surface_reconstruction_boundary_iterator { +private: + Advancing_front_surface_reconstruction_boundary_iterator(){} +public: + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef Advancing_front_surface_reconstruction_boundary_iterator Self; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Vertex Vertex; + + const Surface& S; + int mark; + Finite_vertices_iterator first_vertex; + Vertex_handle pos; + bool first, last; + + Advancing_front_surface_reconstruction_boundary_iterator(const Surface& S_, int m) + : S(S_), mark(m), first_vertex(S.triangulation().finite_vertices_begin()), pos(first_vertex) + { + if (pos->number_of_incident_border() <= 0){ + advance_to_next_boundary(); + } + first = true; + last = false; + } + + Advancing_front_surface_reconstruction_boundary_iterator(const Surface& S_) + : S(S_), pos(NULL) + {} + + Advancing_front_surface_reconstruction_boundary_iterator(const Self& s) + : S(s.S), mark(s.mark), first_vertex(s.first_vertex), pos(s.pos), first(s.first), last(s.last) + {} + + bool operator==(const Self &s) const + { + return pos == s.pos; + } + + bool operator!=(const Self &s) const + { + return pos != s.pos; + } + + + Self operator++() + { + if(pos == NULL) { + return *this; + } + if(first){ + advance_on_boundary(); + first = false; + } else if (last) { + advance_to_next_boundary(); + first = true; + last = false; + } else { + advance_on_boundary(); + if(&*pos == &*first_vertex){ + last = true; + } + } + return *this; + } + + Vertex_handle operator*() + { + return pos; + } + + void advance_on_boundary() + { + if(pos == NULL) { + return; + } + pos = pos->first_incident()->first; + pos->set_post_mark(mark); + } + + void advance_to_next_boundary() + { + if(pos == NULL) { + return; + } + do { + first_vertex++; + } while((first_vertex != S.triangulation().finite_vertices_end()) && + (! ((first_vertex->is_on_border()) + && ! first_vertex->is_post_marked(mark)))); + if(first_vertex != S.triangulation().finite_vertices_end()) { + pos = first_vertex; + pos->set_post_mark(mark); + + } else { + pos = NULL; + } + } +}; + + +template +class Advancing_front_surface_reconstruction { + +public: + typedef Triangulation Triangulation_3; + typedef Advancing_front_surface_reconstruction Extract; + typedef typename Triangulation_3::Geom_traits Geom_traits; + + typedef typename Kernel::FT coord_type; + + typedef typename Kernel::Point_3 Point; + typedef typename Kernel::Vector_3 Vector; + typedef typename Kernel::Segment_3 Segment; + typedef typename Kernel::Triangle_3 Triangle; + typedef typename Kernel::Sphere_3 Sphere; + + typedef typename Triangulation_3::Cell Cell; + typedef typename Triangulation_3::Vertex Vertex; + typedef typename Triangulation_3::Edge Edge; + typedef typename Triangulation_3::Facet Facet; + typedef typename Triangulation_3::Cell_handle Cell_handle; + typedef typename Triangulation_3::Vertex_handle Vertex_handle; + + typedef typename Triangulation_3::Cell_circulator Cell_circulator; + typedef typename Triangulation_3::Facet_circulator Facet_circulator; + + typedef typename Triangulation_3::Locate_type Locate_type; + + typedef typename Triangulation_3::Finite_cells_iterator Finite_cells_iterator; + typedef typename Triangulation_3::Finite_facets_iterator Finite_facets_iterator; + typedef typename Triangulation_3::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Triangulation_3::Finite_edges_iterator Finite_edges_iterator; + + typedef typename Triangulation_3::All_cells_iterator All_cells_iterator; + typedef typename Triangulation_3::All_facets_iterator All_facets_iterator; + typedef typename Triangulation_3::All_vertices_iterator All_vertices_iterator; + typedef typename Triangulation_3::All_edges_iterator All_edges_iterator; + + typedef typename Triangulation_3::Vertex::Edge_incident_facet Edge_incident_facet; + typedef typename Triangulation_3::Vertex::IO_edge_type IO_edge_type; + typedef typename Triangulation_3::Vertex::criteria criteria; + typedef typename Triangulation_3::Vertex::Radius_edge_type Radius_edge_type; + typedef typename Triangulation_3::Vertex::Border_elt Border_elt; + typedef typename Triangulation_3::Vertex::Next_border_elt Next_border_elt; + typedef typename Triangulation_3::Vertex::Radius_ptr_type Radius_ptr_type; + + typedef typename Triangulation_3::Vertex::Incidence_request_iterator Incidence_request_iterator; + typedef typename Triangulation_3::Vertex::Incidence_request_elt Incidence_request_elt; + + typedef std::pair< Vertex_handle, Vertex_handle > Edge_like; + typedef CGAL::Triple< Vertex_handle, Vertex_handle, Vertex_handle > Facet_like; + + + typedef std::multimap< criteria, IO_edge_type*, + std::less > Ordered_border_type; + typedef typename Ordered_border_type::iterator Ordered_border_iterator; + + enum Validation_case {NOT_VALID, NOT_VALID_CONNECTING_CASE, FINAL_CASE, + EAR_CASE, EXTERIOR_CASE, CONNECTING_CASE}; + + //===================================================================== + //===================================================================== +private: + + Triangulation_3& T; + + Ordered_border_type _ordered_border; + int _number_of_border; + + const coord_type SLIVER_ANGULUS; // = sampling quality of the surface + coord_type DELTA; // = sampling quality of the border + coord_type K, min_K; + const coord_type eps; + const coord_type inv_eps_2; // 1/(eps^2) + const coord_type eps_3; // test de ^3 donc points tel 1e-7 soit petit + const criteria STANDBY_CANDIDATE; + const criteria STANDBY_CANDIDATE_BIS; + const criteria NOT_VALID_CANDIDATE; + + const double area; + const double perimeter; + const double abs_area; + const double abs_perimeter; + double total_area; + double total_perimeter; + //--------------------------------------------------------------------- + + CGAL::Timer t1; + + //--------------------------------------------------------------------- + //Pour une visu correcte + //pour retenir les facettes selectionnees + int _vh_number; + int _facet_number; + + //--------------------------------------------------------------------- + //Pour le post traitement + mutable int _postprocessing_counter; + int _size_before_postprocessing; + + std::list outliers; + + int _number_of_connected_components; + +public: + Advancing_front_surface_reconstruction(Triangulation_3& T_, const AFSR_options& opt) + : T(T_), _number_of_border(1), SLIVER_ANGULUS(.86), DELTA(opt.delta), min_K(HUGE_VAL), + eps(1e-7), inv_eps_2(coord_type(1)/(eps*eps)), eps_3(eps*eps*eps), + STANDBY_CANDIDATE(3), STANDBY_CANDIDATE_BIS(STANDBY_CANDIDATE+1), + NOT_VALID_CANDIDATE(STANDBY_CANDIDATE+2), area(opt.area), perimeter(opt.perimeter), + abs_area(opt.abs_area), abs_perimeter(opt.abs_perimeter), total_area(0), total_perimeter(0), + _vh_number(static_cast(T.number_of_vertices())), _facet_number(0), + _postprocessing_counter(0), _size_before_postprocessing(0), _number_of_connected_components(0) + { + bool re_init = false; + do + { + _number_of_connected_components++; + + if ( (re_init = init(re_init)) ) + { + std::cerr << "Growing connected component " << _number_of_connected_components << std::endl; + extend(opt.K_init, opt.K_step, opt.K); + + if ((number_of_facets() > static_cast(T.number_of_vertices()))&& + (opt.NB_BORDER_MAX > 0)) + // en principe 2*nb_sommets = nb_facettes: y a encore de la marge!!! + { + while(postprocessing(opt.NB_BORDER_MAX)){ + extend(opt.K_init, opt.K_step, opt.K); + } + } + } + }while(re_init && + ((_number_of_connected_components < opt.max_connected_comp)|| + (opt.max_connected_comp < 0))); + + _tds_2_inf = AFSR::construct_surface(_tds_2, *this); + } + + ~Advancing_front_surface_reconstruction() + {} + + + typedef Triangulation_data_structure_2, + AFSR::Surface_face_base_2 > TDS_2; + + + mutable TDS_2 _tds_2; + + mutable typename TDS_2::Vertex_handle _tds_2_inf; + + const TDS_2& tds_2() const + { + return _tds_2; + } + + bool is_on_surface(const typename TDS_2::Vertex_handle vh) const + { + return vh != _tds_2_inf; + } + + + + int number_of_connected_components() const + { + return _number_of_connected_components; + } + + + Triangulation_3& + triangulation() const + { + return T; + } + + + + + + int number_of_facets() const + { + return _facet_number; + } + + int number_of_vertices() const + { + return _vh_number; + } + + int number_of_outliers() const + { + return static_cast(outliers.size()); + } + + typedef typename std::list::const_iterator Outlier_iterator; + + Outlier_iterator outliers_begin() const + { + return outliers.begin(); + } + + Outlier_iterator outliers_end() const + { + return outliers.end(); + } + + typedef Advancing_front_surface_reconstruction_boundary_iterator Boundary_iterator; + + Boundary_iterator boundaries_begin() const + { + return Boundary_iterator(*this, get_next_mark()); + } + + Boundary_iterator boundaries_end() const + { + return Boundary_iterator(*this); + } + + + int get_next_mark() const + { + _postprocessing_counter++; + return _postprocessing_counter; + } + + + + Next_border_elt* get_border_elt(const Vertex_handle& v1, const Vertex_handle& v2) + { + return v1->get_border_elt(v2); + } + + //public + + IO_edge_type* get_border_IO_elt(const Vertex_handle& v1, const Vertex_handle& v2) + { + return &get_border_elt(v1,v2)->second.first.second; + } + + IO_edge_type* set_border_elt(const Vertex_handle& v1, const Vertex_handle& v2, + const Border_elt& e) + { + v1->set_next_border_elt(Next_border_elt (v2, e)); + return get_border_IO_elt(v1, v2); + } + + + IO_edge_type* set_again_border_elt(const Vertex_handle& v1, const Vertex_handle& v2, + const Border_elt& e) + { + get_border_elt(v1,v2)->second = e; + return get_border_IO_elt(v1, v2); + } + + //--------------------------------------------------------------------- + + bool is_border_elt(Edge_like& key, Border_elt& result) const + { + Next_border_elt* it12 = key.first->get_border_elt(key.second); + if (it12 != NULL) + { + result = it12->second; + return true; + } + + Next_border_elt* it21 = key.second->get_border_elt(key.first); + if (it21 != NULL) + { + result = it21->second; + std::swap(key.first, key.second); + return true; + } + return false; + } + + //--------------------------------------------------------------------- + bool is_border_elt(Edge_like& key) const { + Next_border_elt* it12 = key.first->get_border_elt(key.second); + if (it12 != NULL) + { + return true; + } + + Next_border_elt* it21 = key.second->get_border_elt(key.first); + if (it21 != NULL) + { + std::swap(key.first, key.second); + return true; + } + return false; + } + //--------------------------------------------------------------------- + + bool is_ordered_border_elt(const Edge_like& key, Border_elt& result) const + { + Next_border_elt* it12 = key.first->get_border_elt(key.second); + if (it12 != NULL) + { + result = it12->second; + return true; + } + return false; + } + + //--------------------------------------------------------------------- + + void + remove_border_elt(const Edge_like& ordered_key) + { + ordered_key.first->remove_border_edge(ordered_key.second); + } + + //--------------------------------------------------------------------- + + bool is_ordered_border_elt(const Edge_like& e, + IO_edge_type* &ptr) const + { + Vertex_handle v1 = e.first; + + Next_border_elt* it12 = v1->get_border_elt(e.second); + if (it12 != NULL) + { + ptr = &it12->second.first.second; + return true; + } + return false; + } + + //--------------------------------------------------------------------- + void set_incidence_request(const Vertex_handle& v, + const criteria& value, + const Edge_like& e) + { + Incidence_request_elt incident_elt(value, e); + v->set_incidence_request(incident_elt); + } + + //--------------------------------------------------------------------- + + bool is_interior_edge(const Edge_like& key) const + // pour gerer certaines aretes interieures: a savoir celle encore connectee au + // bord (en fait seule, les aretes interieures reliant 2 bords nous + // interressent...) + { + return (key.first->is_interior_edge(key.second)|| + key.second->is_interior_edge(key.first)); + } + + //--------------------------------------------------------------------- + +#ifdef AFSR_LAZY + + coord_type get_lazy_squared_radius(const Cell_handle& c) + { + if (c->get_lazy_squared_radius() != NULL) + return *(c->get_lazy_squared_radius()); + + c->set_lazy_squared_radius + (squared_radius(c->vertex(0)->point(), + c->vertex(1)->point(), + c->vertex(2)->point(), + c->vertex(3)->point())); + return *(c->get_lazy_squared_radius()); + } + + Point get_lazy_circumcenter(const Cell_handle& c) + { + if (c->get_lazy_circumcenter() != NULL) + return *(c->get_lazy_circumcenter()); + + c->set_lazy_circumcenter + (circumcenter(c->vertex(0)->point(), + c->vertex(1)->point(), + c->vertex(2)->point(), + c->vertex(3)->point())); + return *(c->get_lazy_circumcenter()); + } + +#endif //NOLAZY + + //--------------------------------------------------------------------- + + Edge_incident_facet next(const Edge_incident_facet& e) const + { + Cell_handle c = e.first.first; + int i = e.second; + int i1 = e.first.second, i2 = e.first.third; + int i3 = (6 - e.second - i1 - i2); + + Cell_handle n = c->neighbor(i); + int j1 = n->index(c->vertex(i1)), j2 = n->index(c->vertex(i2)); + int j = n->index(c->vertex(i3)); + return Edge_incident_facet(Edge(n, j1, j2), j); + } + + //--------------------------------------------------------------------- + + Edge_incident_facet previous(const Edge_incident_facet& e) const + { + Cell_handle c = e.first.first; + int i = e.second; + int i1 = e.first.second, i2 = e.first.third; + int i3 = (6 - e.second - i1 - i2); + + Cell_handle n = c->neighbor(i3); + int j1 = n->index(c->vertex(i1)), j2 = n->index(c->vertex(i2)); + int j = n->index(c->vertex(i)); + return Edge_incident_facet(Edge(n, j1, j2), j); + } + + //--------------------------------------------------------------------- + + bool + my_coplanar(const Point& p, const Point& q, + const Point& r, const Point& s) const + { + coord_type qpx = q.x()-p.x(); + coord_type qpy = q.y()-p.y(); + coord_type qpz = q.z()-p.z(); + coord_type rpx = r.x()-p.x(); + coord_type rpy = r.y()-p.y(); + coord_type rpz = r.z()-p.z(); + coord_type spx = s.x()-p.x(); + coord_type spy = s.y()-p.y(); + coord_type spz = s.z()-p.z(); + + coord_type den = CGAL::determinant(qpx,qpy,qpz, + rpx,rpy,rpz, + spx,spy,spz); + return (CGAL_NTS abs(den) < eps_3); + } + + //--------------------------------------------------------------------- + + + bool + my_collinear(const Point& p, const Point& q, const Point& s) const + { + coord_type psx = p.x()-s.x(); + coord_type psy = p.y()-s.y(); + coord_type psz = p.z()-s.z(); + coord_type qsx = q.x()-s.x(); + coord_type qsy = q.y()-s.y(); + coord_type qsz = q.z()-s.z(); + coord_type rsx = psy*qsz-psz*qsy; + coord_type rsy = psz*qsx-psx*qsz; + coord_type rsz = psx*qsy-psy*qsx; + + coord_type den = CGAL::determinant(psx,psy,psz, + qsx,qsy,qsz, + rsx,rsy,rsz); + + return (CGAL_NTS abs(den) < eps_3); + } + + //--------------------------------------------------------------------- + + void + select_facet(const Cell_handle& c, const int& i) + { + if(area!=0){ + total_area += compute_area(c->vertex((i+1)&3)->point(), + c->vertex((i+2)&3)->point(), + c->vertex((i+3)&3)->point()); + } + if(perimeter != 0){ + total_perimeter += compute_perimeter(c->vertex((i+1)&3)->point(), + c->vertex((i+2)&3)->point(), + c->vertex((i+3)&3)->point()); + } + c->select_facet(i); + _facet_number++; + c->set_facet_number(i, _facet_number); + } + + + void + unselect_facet(const Cell_handle& c, const int& i) + { + c->unselect_facet(i); + } + + int number_of_border_edges() + { + int _border_count(0); + for(Finite_edges_iterator e_it=T.finite_edges_begin(); + e_it!=T.finite_edges_end(); + e_it++) + { + Cell_handle c = (*e_it).first; + int i1 = (*e_it).second, i2 = (*e_it).third; + Edge_like key(c->vertex(i1), c->vertex(i2)); + + if (is_border_elt(key)) + _border_count++; + } + return _border_count; + } + + + + + //===================================================================== + + coord_type get_smallest_radius_delaunay_sphere(const Cell_handle& c, + const int& index) const + { + int i1, i2, i3; + + Cell_handle n = c->neighbor(index); + // lazy evaluation ... + coord_type value = c->get_smallest_radius(index); + if ((value >= 0)&&(n->get_smallest_radius(n->index(c)) == value)) + return value; + + const Point& cp0 = c->vertex(index)->point(); + const Point& cp1 = c->vertex((index+1) & 3)->point(); + const Point& cp2 = c->vertex((index+2) & 3)->point(); + const Point& cp3 = c->vertex((index+3) & 3)->point(); + + const Point& np0 = n->vertex(0)->point(); + const Point& np1 = n->vertex(1)->point(); + const Point& np2 = n->vertex(2)->point(); + const Point& np3 = n->vertex(3)->point(); + + bool c_is_plane(my_coplanar(cp0, cp1, cp2, cp3)); + bool n_is_plane(my_coplanar(np0, np1, np2, np3)); + + bool c_is_infinite(T.is_infinite(c)); + bool n_is_infinite(T.is_infinite(n)); + if ((c_is_plane && n_is_plane)|| + (c_is_plane && n_is_infinite)|| + (n_is_plane && c_is_infinite)|| + my_collinear(cp1, cp2, cp3)) + value = HUGE_VAL; + else + { + if (c_is_infinite||n_is_infinite||c_is_plane||n_is_plane) + { + int ind; + Cell_handle cc; + if(c_is_infinite||c_is_plane) + { + cc = n; + ind = n->index(c); + } + else + { + cc = c; + ind = index; + } + i1 = (ind+1) & 3; + i2 = (ind+2) & 3; + i3 = (ind+3) & 3; + + const Point& pp0 = cc->vertex(ind)->point(); + const Point& pp1 = cc->vertex(i1)->point(); + const Point& pp2 = cc->vertex(i2)->point(); + const Point& pp3 = cc->vertex(i3)->point(); + + Sphere facet_sphere(pp1, pp2, pp3); + if (squared_distance(facet_sphere.center(), pp0) < + facet_sphere.squared_radius()) + { +#ifdef AFSR_LAZY + value = get_lazy_squared_radius(cc); +#else + value = squared_radius(pp0, pp1, pp2, pp3); +#endif + } + else + value = facet_sphere.squared_radius(); + } + else + { + Point cc, cn; +#ifdef AFSR_LAZY + cc = get_lazy_circumcenter(c); + cn = get_lazy_circumcenter(n); +#else + cc = circumcenter(cp0, cp1, cp2, cp3); + cn = circumcenter(np0, np1, np2, np3); +#endif + // computation of the distance of cp1 to the dual segment cc, cn... + Vector V(cc - cn), Vc(cc - cp1), Vn(cp1 - cn); + coord_type ac(V * Vc), an(V * Vn), norm_V(V * V); + if ((ac > 0) && (an > 0)) + { + value = (Vc*Vc) - ac*ac/norm_V; + if ((value < 0)||(norm_V > inv_eps_2)) + value = squared_radius(cp1, cp2, cp3); + } + else + { + if (ac <= 0) + value = squared_distance(cc, cp1); + else // (an <= 0) + value = squared_distance(cn, cp1); + } + } + } + // cache computed values + c->set_smallest_radius(index, value); + n->set_smallest_radius(n->index(c), value); + + return value; + } +//--------------------------------------------------------------------- + coord_type + compute_area(const Point& p1, const Point& p2, const Point& p3) const + { + return typename Kernel::Compute_area_3()(p1,p2,p3); + } +//--------------------------------------------------------------------- + coord_type + compute_perimeter(const Point& p1, const Point& p2, const Point& p3) const + { + return sqrt(squared_distance(p1,p2)) + + sqrt(squared_distance(p2,p3)) + + sqrt(squared_distance(p3,p1)); + } + + //--------------------------------------------------------------------- + // For a border edge e we determine the incident facet which has the highest + // chance to be a natural extension of the surface + + Radius_edge_type compute_value(const Edge_incident_facet& e) + { + Cell_handle c = e.first.first; + int i = e.second; + int i1 = e.first.second, i2 = e.first.third; + int i3 = 6 - e.second - i1 - i2; + + Edge_incident_facet e_it = e, predone = previous(e); + Cell_handle c_predone = predone.first.first; + + coord_type min_valueP = NOT_VALID_CANDIDATE, min_valueA = HUGE_VAL; + Facet min_facet, min_facetA; + bool border_facet(false); + + coord_type pscal; + const Point& p1 = c->vertex(i1)->point(); + const Point& p2 = c->vertex(i2)->point(); + const Point& pc = c->vertex(i3)->point(); + + Vector P2P1 = p1-p2, P2Pn, PnP1; + + Vector v2, v1 = cross_product(pc-p2, P2P1); + + coord_type norm, norm1 = v1*v1; + coord_type norm12 = P2P1*P2P1; + + e_it = next(e_it); + bool succ_start(true); + + do + { + Cell_handle neigh = e_it.first.first; + Facet facet_it(neigh, e_it.second); + + if (!T.is_infinite(facet_it)) + { + int n_ind = facet_it.second; + int n_i1 = e_it.first.second; + int n_i2 = e_it.first.third; + int n_i3 = 6 - n_ind - n_i1 - n_i2; + + coord_type tmp=0; + coord_type ar=0, pe=0; + + // If the triangle has a very large area and a very high perimeter, + // we do not want to consider it as a good candidate. + + if( (abs_area != 0) || ( (_facet_number > 1000) && (area != 0)) ){ + ar = compute_area(facet_it.first->vertex(n_i1)->point(), + facet_it.first->vertex(n_i2)->point(), + facet_it.first->vertex(n_i3)->point()); + } + + if( (abs_perimeter != 0) || ( (_facet_number > 1000) && (perimeter != 0)) ){ + pe = compute_perimeter(facet_it.first->vertex(n_i1)->point(), + facet_it.first->vertex(n_i2)->point(), + facet_it.first->vertex(n_i3)->point()); + } + + if( ((abs_area != 0)&&(ar > abs_area)) || ((abs_perimeter != 0) && (pe > abs_perimeter)) ){ + tmp = HUGE_VAL; + } + + if((tmp != HUGE_VAL) && (_facet_number > 1000)){ + if(area != 0){ + coord_type avg_area = total_area/_facet_number; + if(ar > (area * avg_area)){ + tmp = HUGE_VAL; + } + } + if((perimeter != 0) && (tmp != HUGE_VAL)){ + coord_type avg_perimeter = total_perimeter/_facet_number; + if(pe > (perimeter * avg_perimeter)){ + tmp = HUGE_VAL; + } + } + } + + + if(tmp != HUGE_VAL){ + tmp = get_smallest_radius_delaunay_sphere(neigh, n_ind); + } + + Edge_like el1(neigh->vertex(n_i1),neigh->vertex(n_i3)), + el2(neigh->vertex(n_i2),neigh->vertex(n_i3)); + + if ((tmp != HUGE_VAL)&& + neigh->vertex(n_i3)->not_interior()&& + (!is_interior_edge(el1))&&(!is_interior_edge(el2))) + { + const Point& pn = neigh->vertex(n_i3)->point(); + + P2Pn = pn-p2; + v2 = cross_product(P2P1,P2Pn); + + //pas necessaire de normer pour un bon echantillon: + // on peut alors tester v1*v2 >= 0 + norm = sqrt(norm1 * (v2*v2)); + pscal = v1*v2; + //SLIVER_ANGULUS represente la qualite d'echantillonnage de la + //surface + bool sliver_facet = ((succ_start || (neigh == c_predone))&& + (pscal <= -SLIVER_ANGULUS*norm)); + + if (succ_start) succ_start = false; + + if (!sliver_facet) + { + if (tmp < min_valueA) + { + PnP1 = p1-pn; + // DELTA represente la qualite d'echantillonnage du bord + border_facet = !((P2P1*P2Pn >= + -DELTA*sqrt(norm12*(P2Pn*P2Pn)))&& + (P2P1*PnP1 >= + -DELTA*sqrt(norm12*(PnP1*PnP1)))); + + min_facetA = facet_it; + min_valueA = tmp; + min_valueP = pscal/norm; + } + } + } + } + e_it = next(e_it); + } + while(e_it.first.first != c); + + criteria value; + + if ((min_valueA == HUGE_VAL) || border_facet) // bad facets case + { + min_facet = Facet(c, i); // !!! sans aucune signification.... + value = NOT_VALID_CANDIDATE; // Attention a ne pas inserer dans PQ + } + else + { + min_facet = min_facetA; + + //si on considere seulement la pliure value appartient a [0, 2] + //value = coord_type(1) - min_valueP; + + // si la pliure est bonne on note suivant le alpha sinon on prend en compte la + // pliure seule... pour discriminer entre les bons slivers... + // si on veut discriminer les facettes de bonnes pliures plus finement + // alors -(1+1/min_valueA) app a [-inf, -1] + // -min_valueP app a [-1, 1] + + if (min_valueP > SLIVER_ANGULUS) + value = -(coord_type(1) + coord_type(1)/min_valueA); + else + { + //on refuse une trop grande non-uniformite + coord_type tmp = get_smallest_radius_delaunay_sphere(c, i); + if (min_valueA <= K * tmp) + value = - min_valueP; + else + { + value = STANDBY_CANDIDATE; // tres mauvais candidat mauvaise pliure + // + grand alpha... a traiter plus tard.... + min_K = + (std::min)(min_K, + min_valueA/tmp); + } + } + } + + Cell_handle n = min_facet.first; + int ni1 = n->index(c->vertex(i1)), ni2 = n->index(c->vertex(i2)); + + return + Radius_edge_type(value, IO_edge_type(e, Edge_incident_facet + (Edge(n, ni1, ni2), + min_facet.second))); + } + + //===================================================================== + //===================================================================== + + // The parameter re_init is true the first time only + // Returns true, iff it found a face where the next surface can grow + bool + init(const bool& re_init) + { + Facet min_facet; + coord_type min_value = HUGE_VAL; + int i1, i2, i3; + + if (!re_init) + for(Finite_facets_iterator facet_it = T.finite_facets_begin(); + facet_it != T.finite_facets_end(); + facet_it++) + { + coord_type value = get_smallest_radius_delaunay_sphere((*facet_it).first, + (*facet_it).second); + if (value < min_value) + { + min_facet = *facet_it; + min_value = value; + } + } + else //if (re_init) + for(Finite_facets_iterator facet_it = T.finite_facets_begin(); + facet_it != T.finite_facets_end(); + facet_it++) + { + Cell_handle c = (*facet_it).first; + int index = (*facet_it).second; + if (c->vertex((index+1) & 3)->is_exterior()) + if (c->vertex((index+2) & 3)->is_exterior()) + if (c->vertex((index+3) & 3)->is_exterior()) + { + coord_type value = get_smallest_radius_delaunay_sphere(c, index); + coord_type pe=0; + if(abs_perimeter != 0){ + pe = compute_perimeter(c->vertex((index+1)&3)->point(), + c->vertex((index+2)&3)->point(), + c->vertex((index+3)&3)->point()); + if(pe > abs_perimeter){ + value = min_value; + } + } + + if (value < min_value) + { + min_facet = *facet_it; + min_value = value; + } + } + } + + if (min_value != HUGE_VAL) + { + Cell_handle c_min = min_facet.first; + + int ind = min_facet.second; + i1 = (ind+1) & 3; + i2 = (ind+2) & 3; + i3 = (ind+3) & 3; + + Radius_edge_type e12, e23, e31; + + e12 = compute_value(Edge_incident_facet(Edge(c_min, i1, i2), ind)); + e23 = compute_value(Edge_incident_facet(Edge(c_min, i2, i3), ind)); + e31 = compute_value(Edge_incident_facet(Edge(c_min, i3, i1), ind)); + + IO_edge_type* p12 = set_border_elt(c_min->vertex(i1), c_min->vertex(i2), + Border_elt(e12, _number_of_border)); + IO_edge_type* p23 = set_border_elt(c_min->vertex(i2), c_min->vertex(i3), + Border_elt(e23, _number_of_border)); + IO_edge_type* p31 = set_border_elt(c_min->vertex(i3), c_min->vertex(i1), + Border_elt(e31, _number_of_border)); + + c_min->vertex(i1)->inc_mark(); + c_min->vertex(i2)->inc_mark(); + c_min->vertex(i3)->inc_mark(); + _ordered_border.insert(Radius_ptr_type (e12.first, p12)); + _ordered_border.insert(Radius_ptr_type (e23.first, p23)); + _ordered_border.insert(Radius_ptr_type (e31.first, p31)); + + select_facet(c_min, ind); + return true; + } + return false; + } + + //--------------------------------------------------------------------- + // test de reciprocite avant de recoller une oreille anti-singularite + int + test_merge(const Edge_like& ordered_key, const Border_elt& result, + const Vertex_handle& v, const coord_type& ear_alpha) + { + Edge_incident_facet Ifacet = result.first.second.first; + + const Point& p1 = (ordered_key.first)->point(); + const Point& p2 = (ordered_key.second)->point(); + const Point& pc = v->point(); + + Cell_handle neigh = Ifacet.first.first; + int n_ind = Ifacet.second; + int n_i1 = Ifacet.first.second; + int n_i2 = Ifacet.first.third; + int n_i3 = (6 - n_ind - n_i1 - n_i2); + + const Point& pn = neigh->vertex(n_i3)->point(); + Vector v1 = cross_product(pc-p2,p1-p2), + v2 = cross_product(p1-p2,pn-p2); + coord_type norm = sqrt((v1*v1)*(v2*v2)); + + if (v1*v2 > SLIVER_ANGULUS*norm) + return 1; // label bonne pliure sinon: + + if (ear_alpha <= K*get_smallest_radius_delaunay_sphere(neigh, n_ind)) + return 2; // label alpha coherent... + + return 0; //sinon oreille a rejeter... + } + + + //--------------------------------------------------------------------- + + void + ordered_map_erase(const criteria& value, const IO_edge_type* pkey) + { + std::size_t number_of_conflict = _ordered_border.count(value); + if (number_of_conflict == 1) + { + _ordered_border.erase(_ordered_border.find(value)); + + } + + else if (number_of_conflict > 1) + { + Ordered_border_iterator elt_it = + _ordered_border.find(value); + // si ca foire jamais on peut s'areter des que l'elt + // est trouve!!! + for(std::size_t jj=0; (jjsecond) == ((long) pkey)) + { + _ordered_border.erase(elt_it); + return; + } + elt_it++; + } + } + } + + //--------------------------------------------------------------------- + + void + force_merge(const Edge_like& ordered_key, const Border_elt& result) + { + criteria value = result.first.first; + IO_edge_type* pkey = get_border_IO_elt(ordered_key.first, ordered_key.second); + + ordered_map_erase(value, pkey); + + remove_border_elt(ordered_key); + } + + //--------------------------------------------------------------------- + + void dequeue_incidence_request(const Vertex_handle& v) + { + if (v->is_incidence_requested()) + { + for(Incidence_request_iterator v_it = v->incidence_request_begin(); + v_it != v->get_incidence_request_end(); + v_it++) + { + IO_edge_type* ptr; + + if (is_ordered_border_elt(v_it->second, ptr)) + _ordered_border.insert(Radius_ptr_type(v_it->first, ptr)); + } + v->erase_incidence_request(); + } + } + + + //--------------------------------------------------------------------- + + void + merge_ear(const Edge_like& ordered_el1, const Border_elt& result1, + const Edge_like& ordered_key, + const Vertex_handle& v1, const Vertex_handle& v2, + const Edge_incident_facet& edge_Ifacet_2) + { + remove_border_elt(ordered_key); + force_merge(ordered_el1, result1); + + Radius_edge_type e2 = compute_value(edge_Ifacet_2); + + IO_edge_type* p2; + if (ordered_el1.first == v1) + p2 = set_border_elt(v2, ordered_el1.second, + Border_elt(e2,result1.second)); + else + p2 = set_border_elt(ordered_el1.first, v2, + Border_elt(e2,result1.second)); + + v1->dec_mark(); + + _ordered_border.insert(Radius_ptr_type(e2.first, p2)); + + //depiler les eventuelles requettes de connections avortees... zones etoilees, + //en effet le bord a change donc on peut peut etre maintenant. + dequeue_incidence_request(v2); + + if (ordered_el1.first == v1) + dequeue_incidence_request(ordered_el1.second); + else + dequeue_incidence_request(ordered_el1.first); + } + + //--------------------------------------------------------------------- + + void + border_extend(const Edge_like& ordered_key, const Border_elt& result12, + const Vertex_handle& v1, const Vertex_handle& v2, + const Vertex_handle& v3, + const Radius_edge_type& e1, const Radius_edge_type& e2, + IO_edge_type* &p1, IO_edge_type* &p2) + { + remove_border_elt(ordered_key); + + //depiler v3 avant de le mettre a jour... pour reperer s'il est sur un bord + if (v3->is_on_border()) + dequeue_incidence_request(v3); + + if (ordered_key.first == v1) + { + p1 = set_border_elt(v1, v3, Border_elt(e1,result12.second)); + p2 = set_border_elt(v3, v2, Border_elt(e2,result12.second)); + } + else + { + p2 = set_border_elt(v2, v3, Border_elt(e2,result12.second)); + p1 = set_border_elt(v3, v1, Border_elt(e1,result12.second)); + } + + v3->inc_mark(); + + //depiler les eventuelles requettes de connections avortees... zones etoilees, + //en effet le bord a change donc on peut peut etre maintenant. + dequeue_incidence_request(v1); + dequeue_incidence_request(v2); + } + + //===================================================================== + + Validation_case + validate(const Edge_incident_facet& edge_Efacet, + const criteria& value) + { + int i = (6 - edge_Efacet.second + - edge_Efacet.first.second + - edge_Efacet.first.third); + Cell_handle c = edge_Efacet.first.first; + + Vertex_handle v1 = c->vertex(edge_Efacet.first.second), + v2 = c->vertex(edge_Efacet.first.third); + + Edge_like ordered_el1(c->vertex(i), v1); + Edge_like ordered_el2(c->vertex(i), v2); + Border_elt result1, result2, result12; + + Edge_like ordered_key(v1,v2); + + if (!is_border_elt(ordered_key, result12)) + std::cerr << "+++probleme coherence bord " << std::endl; + + bool is_border_el1 = is_border_elt(ordered_el1, result1), + is_border_el2 = is_border_elt(ordered_el2, result2); + + Radius_edge_type e1, e2; + + if (c->vertex(i)->not_interior()) + { + if ((!is_interior_edge(ordered_el1))&& + (!is_interior_edge(ordered_el2))) + { + //toujours utile meme avec l'essai de try_to_close_border avant + //validate pour la resolution de singularite par oreille qui elle + //doit etre dans Delaunay. + if (is_border_el1&&is_border_el2) + { + remove_border_elt(ordered_key); + force_merge(ordered_el1, result1); + force_merge(ordered_el2, result2); + + v1->dec_mark(); + v2->dec_mark(); + c->vertex(i)->dec_mark(); + + select_facet(c, edge_Efacet.second); + + return FINAL_CASE; + } + + //--------------------------------------------------------------------- + //on peut alors marquer v1 et on pourrait essayer de merger + //sans faire de calcul inutile??? + if (is_border_el1) + { + Edge_incident_facet edge_Ifacet_2(Edge(c, i, edge_Efacet.first.third), + edge_Efacet.second); + merge_ear(ordered_el1, result1, + ordered_key, v1, v2, edge_Ifacet_2); + + select_facet(c, edge_Efacet.second); + + return EAR_CASE; + } + + //--------------------------------------------------------------------- + //idem pour v2 + if (is_border_el2) + { + Edge_incident_facet edge_Ifacet_1(Edge(c, i, edge_Efacet.first.second), + edge_Efacet.second); + merge_ear(ordered_el2, result2, + ordered_key, v2, v1, edge_Ifacet_1); + + select_facet(c, edge_Efacet.second); + + return EAR_CASE; + } + + + //--------------------------------------------------------------------- + if ((!is_border_el1)&&(!is_border_el2)) + { + // si on veut s'interdir de spliter un bord (pelure d'orange....) + // seulement c->vertex(i)->is_exterior() + // pour s'autoriser des split de bord surface a bord->sphere ou Moebius... + // alors || is_on_same_border: + // if (c->vertex(i)->is_exterior() || is_on_same_border) + // pour passer au tore (changementde type de topologie) + // recoller deux bord different... + // if (c->vertex(i)->not_interior() deja teste en haut + + if(c->vertex(i)->is_exterior()) + { + Edge_incident_facet edge_Ifacet_1(Edge(c, i, edge_Efacet.first.second), + edge_Efacet.second); + + Edge_incident_facet edge_Ifacet_2(Edge(c, i, edge_Efacet.first.third), + edge_Efacet.second); + e1 = compute_value(edge_Ifacet_1); + e2 = compute_value(edge_Ifacet_2); + + IO_edge_type* p1; + IO_edge_type* p2; + + border_extend(ordered_key, result12, + v1, v2, c->vertex(i), + e1, e2, p1, p2); + + // if e1 contain HUGE_VAL there is no candidates to + // continue: compute_value is not valid... + + _ordered_border.insert(Radius_ptr_type(e1.first, p1)); + + _ordered_border.insert(Radius_ptr_type(e2.first, p2)); + + select_facet(c, edge_Efacet.second); + + return EXTERIOR_CASE; + } + else // c->vertex(i) is a border point (and now there's only 1 + // border incident to a point... _mark<1 even if th orientation + // may be such as one vh has 2 successorson the same border... + { + // a ce niveau on peut tester si le recollement se fait en + // maintenant la compatibilite d'orientation des bords (pour + // surface orientable...) ou si elle est brisee... + Edge_incident_facet edge_Ifacet_1(Edge(c, i, edge_Efacet.first.second), + edge_Efacet.second); + Edge_incident_facet edge_Ifacet_2(Edge(c, i, edge_Efacet.first.third), + edge_Efacet.second); + + e1 = compute_value(edge_Ifacet_1); + e2 = compute_value(edge_Ifacet_2); + + if ((e1.first >= STANDBY_CANDIDATE)&&(e2.first >= STANDBY_CANDIDATE)) + return NOT_VALID_CONNECTING_CASE; + + // vu compute value: les candidats oreilles fournis sont sans + // aretes interieures et le sommet oppose n'est pas non plus interieur + Edge_incident_facet ear1 = e1.second.second; + Edge_incident_facet ear2 = e2.second.second; + + int ear1_i = (6 - ear1.second + - ear1.first.second + - ear1.first.third); + Cell_handle ear1_c = ear1.first.first; + Border_elt result_ear1; + + int ear2_i = (6 - ear2.second + - ear2.first.second + - ear2.first.third); + Cell_handle ear2_c = ear2.first.first; + Border_elt result_ear2; + + Edge_like ear1_e, ear2_e; + // pour maintenir la reconstruction d'une surface orientable : + // on verifie que les bords se recollent dans des sens opposes + if (ordered_key.first==v1) + { + ear1_e = Edge_like(c->vertex(i), ear1_c ->vertex(ear1_i)); + ear2_e = Edge_like(ear2_c ->vertex(ear2_i), c->vertex(i)); + } + else + { + ear1_e = Edge_like(ear1_c ->vertex(ear1_i), c->vertex(i)); + ear2_e = Edge_like(c->vertex(i), ear2_c ->vertex(ear2_i)); + } + + //maintient la surface orientable + bool is_border_ear1 = is_ordered_border_elt(ear1_e, result_ear1); + bool is_border_ear2 = is_ordered_border_elt(ear2_e, result_ear2); + bool ear1_valid(false), ear2_valid(false); + + if (is_border_ear1&&(e1.first < STANDBY_CANDIDATE)&& + (e1.first <= value)&& + (result12.second==result_ear1.second)) + { + ear1_valid = test_merge(ear1_e, result_ear1, v1, + get_smallest_radius_delaunay_sphere(ear1_c, + ear1.second)) != 0; + } + + if (is_border_ear2&&(e2.first < STANDBY_CANDIDATE)&& + (e2.first <= value)&& + (result12.second==result_ear2.second)) + { + ear2_valid = test_merge(ear2_e, result_ear2, v2, + get_smallest_radius_delaunay_sphere(ear2_c, + ear2.second)) != 0; + } + + if ((!ear1_valid)&&(!ear2_valid)) + return NOT_VALID_CONNECTING_CASE; + + IO_edge_type* p1; + IO_edge_type* p2; + + border_extend(ordered_key, result12, + v1, v2, c->vertex(i), + e1, e2, p1, p2); + + if (ear1_valid&&ear2_valid&&(ear1_e==ear2_e)) + { + if (e1.first < e2.first) + { + Validation_case res = validate(ear1, e1.first); + if (!((res == EAR_CASE)||(res == FINAL_CASE))) + std::cerr << "+++probleme de recollement : cas " + << res << std::endl; + e2 = compute_value(edge_Ifacet_2); + + if (ordered_key.first == v1) + p2 = set_again_border_elt(c->vertex(i), v2, + Border_elt(e2, result2.second)); + else + p2 = set_again_border_elt(v2, c->vertex(i), + Border_elt(e2, result2.second)); + + _ordered_border.insert(Radius_ptr_type(e2.first, p2)); + } + else + { + Validation_case res = validate(ear2, e2.first); + if (!((res == EAR_CASE)||(res == FINAL_CASE))) + std::cerr << "+++probleme de recollement : cas " + << res << std::endl; + e1 = compute_value(edge_Ifacet_1); + + if (ordered_key.first == v1) + p1 = set_again_border_elt(v1, c->vertex(i), + Border_elt(e1, result1.second)); + else + p1 = set_again_border_elt(c->vertex(i), v1, + Border_elt(e1, result1.second)); + + _ordered_border.insert(Radius_ptr_type(e1.first, p1)); + } + } + else// les deux oreilles ne se recollent pas sur la meme arete... + { + // on resoud la singularite. + if (ear1_valid) + { + Validation_case res = validate(ear1, e1.first); + if (!((res == EAR_CASE)||(res == FINAL_CASE))) + std::cerr << "+++probleme de recollement : cas " + << res << std::endl; + } + if (ear2_valid) + { + Validation_case res = validate(ear2, e2.first); + if (!((res == EAR_CASE)||(res == FINAL_CASE))) + std::cerr << "+++probleme de recollement : cas " + << res << std::endl; + } + // on met a jour la PQ s'il y a lieu... mais surtout pas + // avant la resolution de la singularite + if (!ear1_valid) + { + _ordered_border.insert(Radius_ptr_type(e1.first, p1)); + } + if (!ear2_valid) + { + _ordered_border.insert(Radius_ptr_type(e2.first, p2)); + } + } + + select_facet(c, edge_Efacet.second); + + return CONNECTING_CASE; + } + } + } + } + return NOT_VALID; + } + + //===================================================================== + void re_compute_values() + { + if(!_ordered_border.empty()) + { + Ordered_border_type _ordered_border_tmp; + do + { + Ordered_border_iterator e_it = _ordered_border.begin(); + Edge_incident_facet mem_Ifacet = e_it->second->first; + Cell_handle c_tmp = mem_Ifacet.first.first; + _ordered_border.erase(e_it); + Vertex_handle v1 = c_tmp->vertex(mem_Ifacet.first.second); + Vertex_handle v2 = c_tmp->vertex(mem_Ifacet.first.third); + + Radius_edge_type new_candidate; + new_candidate = compute_value(mem_Ifacet); + + if (new_candidate.first == STANDBY_CANDIDATE) + { + // a garder pour un K un peu plus grand... + new_candidate.first = STANDBY_CANDIDATE_BIS; + } + + Border_elt result; + Edge_like key_tmp(v1,v2); + is_border_elt(key_tmp, result); + IO_edge_type* pnew = + set_again_border_elt(key_tmp.first, key_tmp.second, + Border_elt (new_candidate, result.second)); + _ordered_border_tmp.insert(Radius_ptr_type(new_candidate.first, pnew)); + } + while(!_ordered_border.empty()); + + _ordered_border.swap(_ordered_border_tmp); + } + } + + //--------------------------------------------------------------------- + + void + extend(const coord_type& K_init, const coord_type& K_step, const coord_type& K_max) + { + // initilisation de la variable globale K: qualite d'echantillonnage requise + K = K_init; // valeur d'initialisation de K pour commencer prudemment... + + Vertex_handle v1, v2; + t1.start(); + if (_ordered_border.empty()) return; + do + { + min_K = HUGE_VAL; // pour retenir le prochain K necessaire pour progresser... + do + { + Ordered_border_iterator e_it = _ordered_border.begin(); + + criteria value = e_it->first; + if (value >= STANDBY_CANDIDATE) + re_compute_values(); + else + { + Edge_incident_facet candidate = e_it->second->second; + Cell_handle c_ext = candidate.first.first; + int i1, i2 , i3; + i1 = candidate.first.second; + i2 = candidate.first.third; + i3 = (6 - i1- i2 - candidate.second); + + Edge_incident_facet mem_Ifacet = e_it->second->first; + Cell_handle c_tmp = mem_Ifacet.first.first; + + v1 = c_tmp->vertex(mem_Ifacet.first.second); + v2 = c_tmp->vertex(mem_Ifacet.first.third); + + Radius_edge_type mem_e_it(e_it->first, *e_it->second); + + _ordered_border.erase(e_it); + + Validation_case validate_result = validate(candidate, value); + + if ((validate_result == NOT_VALID)|| + (validate_result == NOT_VALID_CONNECTING_CASE)) + { + Radius_edge_type new_candidate; + Border_elt result; + Edge_like key_tmp(v1,v2); + is_border_elt(key_tmp, result); + + if (validate_result == NOT_VALID_CONNECTING_CASE) + set_incidence_request(c_ext->vertex(i3), value, key_tmp); + + if (validate_result == NOT_VALID) + { + new_candidate = compute_value(mem_Ifacet); + if ((new_candidate != mem_e_it)) + // &&(new_candidate.first < NOT_VALID_CANDIDATE)) + { + IO_edge_type* pnew = + set_again_border_elt(key_tmp.first, key_tmp.second, + Border_elt (new_candidate, result.second)); + _ordered_border.insert(Radius_ptr_type(new_candidate.first, + pnew)); + } + } + } + } + } + while((!_ordered_border.empty())&& + (_ordered_border.begin()->first < STANDBY_CANDIDATE_BIS)); + + K += (std::max)(K_step, min_K-K+eps); + // on augmente progressivement le K mais on a deja rempli sans + // faire des betises auparavant... + } + while((!_ordered_border.empty())&&(K <= K_max)&&(min_K != HUGE_VAL)); + t1.stop(); + +#ifdef VERBOSE + if ((min_K < HUGE_VAL)&&(!_ordered_border.empty())) { + std::cout << " [ next K required = " << min_K << " ]" << std::endl; + } +#endif // VERBOSE + } + + + + + //--------------------------------------------------------------------- + // En principe, si l'allocateur de cellules etait bien fait on aurait pas besoin + // de mettre a jour les valeurs rajoutees pour les cellules a la main... + + void re_init_for_free_cells_cache(const Vertex_handle& vh) + { + std::list ch_set; + T.incident_cells(vh, std::back_inserter(ch_set)); + for (typename std::list::iterator c_it = ch_set.begin(); + c_it != ch_set.end(); + c_it++) + (*c_it)->clear(); + } + + + //--------------------------------------------------------------------- + + void swap_selected_facets_on_conflict_boundary(const Vertex_handle& vh) + { + std::list ch_set; + T.incident_cells(vh, std::back_inserter(ch_set)); + for (typename std::list::iterator c_it = ch_set.begin(); + c_it != ch_set.end(); c_it++) + { + Cell_handle c = *c_it; + int index = c->index(vh); + Cell_handle neigh = c->neighbor(index); + int n_ind = neigh->index(c); + neigh->set_smallest_radius(n_ind, -1); // pour obliger le recalcul + // si c est selectionnee c'est qu'elle est aussi le mem_IFacet renvoye par + // compute_value... donc a swapper aussi + if (c->is_selected_facet(index)) + { + int fn = c->facet_number(index); + unselect_facet(c, index); + neigh->select_facet(n_ind); + neigh->set_facet_number(n_ind, fn); + int i1 = (n_ind+1) & 3; + int i2 = (n_ind+2) & 3; + int i3 = (n_ind+3) & 3; + Edge_like key(neigh->vertex(i1), neigh->vertex(i2)); + + if (is_border_elt(key)) + { + Edge_incident_facet ei_facet(Edge(neigh, i1, i2), + n_ind); + *get_border_IO_elt(key.first, key.second) = + IO_edge_type(ei_facet, ei_facet); + } + key = Edge_like(neigh->vertex(i1), neigh->vertex(i3)); + if (is_border_elt(key)) + { + Edge_incident_facet ei_facet(Edge(neigh, i1, i3), + n_ind); + *get_border_IO_elt(key.first, key.second) = + IO_edge_type(ei_facet, ei_facet); + } + key = Edge_like(neigh->vertex(i3), neigh->vertex(i2)); + if (is_border_elt(key)) + { + Edge_incident_facet ei_facet(Edge(neigh, i3, i2), + n_ind); + *get_border_IO_elt(key.first, key.second) = + IO_edge_type(ei_facet, ei_facet); + } + } + } + } + + //--------------------------------------------------------------------- + + Facet next_surface_facet(const Edge_incident_facet& start) + { + Edge_incident_facet circ = next(start); + Cell_handle c = start.first.first; + do + { + Cell_handle ch = circ.first.first; + int ind = circ.second; + Cell_handle neigh = ch->neighbor(ind); + int n_ind = neigh->index(ch); + if (ch->is_selected_facet(ind)) + return Facet(ch, ind); + if (neigh->is_selected_facet(n_ind)) + return Facet(neigh, n_ind); + circ = next(circ); + } + while(circ.first.first != c); + // si on passe par la, alors y a eu un probleme.... + std::cerr << "+++probleme dans la MAJ avant remove..." << std::endl; + + return Facet(c, start.second); + } + + //--------------------------------------------------------------------- + + void retract_border_for_incident_facets(const Vertex_handle& vh) + { + Next_border_elt border_elt = *(vh->first_incident()); + int border_index = border_elt.second.second; + Vertex_handle vh_succ = border_elt.first; + IO_edge_type io_edge = border_elt.second.first.second; + Edge_incident_facet i_facet = io_edge.first; + Cell_handle c = i_facet.first.first; + int i1 = c->index(vh); + int i2 = c->index(vh_succ); + int index = i_facet.second; + int i3 = 6 - index - i1 - i2; + Vertex_handle vh_int = c->vertex(i3); + ordered_map_erase(border_elt.second.first.first, + get_border_IO_elt(vh, vh_succ)); + vh->remove_border_edge(vh_succ); + // 1- a virer au cas ou car vh va etre detruit + vh_succ->remove_interior_edge(vh); + bool while_cond(true); + do + { + _facet_number--; + + assert(c->is_selected_facet(index)); + unselect_facet(c, index); + + Facet f32 = + next_surface_facet(Edge_incident_facet(Edge(c, i3, i2), + index)); + + if (!vh_int->is_on_border()) + { + vh_int->re_init(); + vh_int->inc_mark(); + } + + Edge_incident_facet e32(Edge(f32.first, + f32.first->index(vh_int), + f32.first->index(vh_succ)), f32.second); + Radius_edge_type rad_elt_32(STANDBY_CANDIDATE, IO_edge_type(e32, e32)); + Border_elt result; + if (is_ordered_border_elt(Edge_like(vh_int, vh), result)) + { + ordered_map_erase(result.first.first, get_border_IO_elt(vh_int, vh)); + vh_int->remove_border_edge(vh); + // 1- a virer au cas ou car vh va etre detruit + vh_int->remove_interior_edge(vh); + while_cond = false; + } + // a titre preventif... on essaye de s'assurer de marquer les aretes + // interieures au sens large... + + // 2- a virer a tout pris pour que maintenir le sens de interior edge + vh_int->remove_interior_edge(vh_succ); + vh_succ->remove_interior_edge(vh_int); + + IO_edge_type* p32 = set_border_elt(vh_int, vh_succ, + Border_elt(rad_elt_32, border_index)); + _ordered_border.insert(Radius_ptr_type (STANDBY_CANDIDATE, p32)); + + // incrementation... + if (while_cond) + { + Facet f31 = + next_surface_facet(Edge_incident_facet(Edge(c, i3, i1), + index)); + + c = f31.first; + index = f31.second; + i1 = c->index(vh); + vh_succ = vh_int; + i2 = c->index(vh_int); + i3 = 6 - index - i1 - i2; + vh_int = c->vertex(i3); + } + } + while(while_cond); + } + + //--------------------------------------------------------------------- + + bool create_singularity(const Vertex_handle& vh) + { + // Pour reperer le cas de triangle isole + if (vh->is_on_border()) + { + // vh sommet 0 + Next_border_elt border_elt = *(vh->first_incident()); + Vertex_handle vh_1 = border_elt.first;// sommet 1 + border_elt = *(vh_1->first_incident()); + Vertex_handle vh_2 = border_elt.first;// sommet 2 + border_elt = *(vh_2->first_incident()); + Vertex_handle vh_3 = border_elt.first;// sommet 0 ??? + Cell_handle c; + int i, j, k; + if ((vh_3 == vh)&&(T.is_facet(vh, vh_1, vh_2, c, i ,j ,k))) + { + int l = 6-i-j-k; + Cell_handle neigh = c->neighbor(l); + + if + (c->is_selected_facet(l)||neigh->is_selected_facet(neigh->index(c))) + return true; + } + } + + + // Reperer le cas d'aretes interieures... + std::list vh_list; + T.incident_vertices(vh, std::back_inserter(vh_list)); + + for (typename std::list::iterator v_it = vh_list.begin(); + v_it != vh_list.end(); v_it++) + if ((*v_it)->is_on_border() && is_interior_edge(Edge_like(vh, *v_it))) + return true; + return false; + } + + + //--------------------------------------------------------------------- + + void + store_outlier(const Point& p){ + outliers.push_back(p); + } + + void dec_vh_number() + { + _vh_number--; + } + + struct Remove : public std::unary_function + { + + Extract& E; + Triangulation_3& T; + + Remove(Extract& E_, Triangulation_3& T_) : E(E_), T(T_) {} + + bool operator()(Vertex_handle vh) { + if (vh->is_exterior()) + { + E.swap_selected_facets_on_conflict_boundary(vh); + E.re_init_for_free_cells_cache(vh); + Point p = vh->point(); + T.remove(vh); + E.dec_vh_number(); + E.store_outlier(p); + + return true; + } + else if (vh->is_on_border()&&(!E.create_singularity(vh))) + { + E.swap_selected_facets_on_conflict_boundary(vh); + E.retract_border_for_incident_facets(vh); + E.re_init_for_free_cells_cache(vh); + Point p = vh->point(); + T.remove(vh); + E.dec_vh_number(); + E.store_outlier(p); + + return true; + } + else + { } + return false; + } + }; + + + //--------------------------------------------------------------------- + + bool postprocessing(int NB_BORDER_MAX) + { + _postprocessing_counter++; + + std::list L_v; + + // Pour controler les sommets choisis sur le bord... + + // nombre d'aretes a partir duquel on considere que c'est irrecuperable NB_BORDER_MAX + + int vh_on_border_inserted(0); + for(Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++) + { + v_it->erase_incidence_request(); + if ((v_it->is_on_border())&& + (!v_it->is_post_marked(_postprocessing_counter))) + { + std::list L_v_tmp; + Vertex_handle vprev_it(v_it), done(vprev_it), vh_it; + // Vertex_handle vsucc_it; + int v_count(0); + // collect all vertices on the border + do + { + vh_it = vprev_it->first_incident()->first; + L_v_tmp.push_back(vh_it); + vh_it->set_post_mark(_postprocessing_counter); + vprev_it = vh_it; + v_count++; + } + while((vprev_it != done)&&(v_count < NB_BORDER_MAX)); + // we stopped either because we did a complete tour, or because + // the border was so long that we consider it as too big to close + // e.g., if it is a terrain with only one real border at the exterior + if (v_count < NB_BORDER_MAX) + { + L_v.insert(L_v.begin(), L_v_tmp.begin(), L_v_tmp.end()); + vh_on_border_inserted += v_count; + } + + } + if (v_it->is_exterior()) + L_v.push_back(v_it); + } + + std::size_t itmp, L_v_size_mem; + L_v_size_mem = L_v.size(); + if ((vh_on_border_inserted != 0)&& // pour ne post-traiter que les bords + (L_v.size() < .1 * _size_before_postprocessing)) + { + { + do + { + itmp = L_v.size(); + typename std::list::iterator new_end = + std::remove_if(L_v.begin(), L_v.end(), Remove(*this,T)); + L_v.erase(new_end, L_v.end()); + } + while (!L_v.empty() && (L_v.size() < itmp)); + } +#ifdef VERBOSE + if(L_v.size() > 0){ + std::cout << " " << L_v.size() << " non regular points." << std::endl; + } +#endif // VERBOSE + re_compute_values(); + } + else{ + return false; + } + // we stop if we removed more than 10% of points or after 20 rounds + if ((L_v_size_mem == L_v.size())|| + ((_size_before_postprocessing - T.number_of_vertices()) > + .1 * _size_before_postprocessing)|| + (_postprocessing_counter > 20)){ + return false; + } + + min_K = HUGE_VAL; + // fin-- + // if (_postprocessing_counter < 5) + // return true; + return true; + } + + +}; // class Advancing_front_surface_reconstruction + +} // namespace CGAL + +#endif // CGAL_ADVANCING_FRONT_SURFACE_RECONSTRUCTION_H diff --git a/Advancing_front_surface_reconstruction/include/CGAL/IO/AFSR_vrml.h b/Advancing_front_surface_reconstruction/include/CGAL/IO/AFSR_vrml.h new file mode 100644 index 00000000000..aa399136818 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/IO/AFSR_vrml.h @@ -0,0 +1,68 @@ +#ifndef CGAL_AFSR_VRML_H +#define CGAL_AFSR_VRML_H + +namespace CGAL { + +template < class Vb, class Fb> +void +afsr_vrml_output(const Triangulation_data_structure_2& tds, + std::ostream& os, double r, double g, double b, + typename Triangulation_data_structure_2::Vertex_handle v, bool skip_infinite) +{ + typedef typename Triangulation_data_structure_2::Vertex_handle Vertex_handle; + typedef typename Triangulation_data_structure_2::Vertex_iterator Vertex_iterator; + typedef typename Triangulation_data_structure_2::Face_iterator Face_iterator; + + // ouput to a vrml file style + // Point are assumed to be 3d points with a stream operator << + // if non NULL, v is the vertex to be output first + // if skip_inf is true, the point in the first vertex is not output + // and the faces incident to v are not output + // (it may be for instance the infinite vertex of the terrain) + + os << "#VRML V2.0 utf8" << std::endl; + os << "Shape {\n" + << "appearance Appearance {\n" + << "material Material { diffuseColor " << r << " " << g << " " << b << "}}\n"; + os << "\tgeometry IndexedFaceSet {" << std::endl; + os << "\t\tcoord Coordinate {" << std::endl; + os << "\t\t\tpoint [" << std::endl; + + std::map vmap; + Vertex_iterator vit; + Face_iterator fit; + + int inum = 0; + for( vit= tds.vertices_begin(); vit != tds.vertices_end() ; ++vit) { + if ( v != vit) { + vmap[vit] = inum++; + os << "\t\t\t\t" << *vit << ","<< std::endl; + } + } + + os << "\t\t\t]" << std::endl; + os << "\t\t}" << std::endl; + os << "\t\tsolid FALSE\n" + "\t\tcoordIndex [" << std::endl; + + // faces + for(fit= tds.faces_begin(); fit != tds.faces_end(); ++fit) { + if (!skip_infinite || !fit->has_vertex(v)) { + os << "\t\t\t"; + os << vmap[(*fit).vertex(0)] << ", "; + os << vmap[(*fit).vertex(1)] << ", "; + os << vmap[(*fit).vertex(2)] << ", "; + os << "-1, " << std::endl; + } + } + os << "\t\t]" << std::endl; + os << "\t}" << std::endl; + os << "}" << std::endl; + return; +} + + +} // namespace CGAL + +#endif + diff --git a/Advancing_front_surface_reconstruction/include/CGAL/IO/Advancing_front_surface_reconstruction.h b/Advancing_front_surface_reconstruction/include/CGAL/IO/Advancing_front_surface_reconstruction.h new file mode 100644 index 00000000000..4adca6f951a --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/IO/Advancing_front_surface_reconstruction.h @@ -0,0 +1,1029 @@ +#ifndef CGAL_IO_ADVANCING_FRONT_SURFACE_RECONSTRUCTION_H +#define CGAL_IO_ADVANCING_FRONT_SURFACE_RECONSTRUCTION_H + + +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace CGAL { + + +template +struct Is_not_exterior { + bool operator()(const Vertex& v)const { + return ! v.is_exterior(); + } +}; + +template +void +write_to_file_medit(char* foutput, const Surface& S) +{ + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Vertex Vertex; + typedef typename Surface::Cell_handle Cell_handle; + + Triangulation_3& T = S.triangulation(); + + char foutput_points[100]; + char foutput_faces[100]; + std::strcpy(foutput_points, foutput); + std::strcpy(foutput_faces, foutput); + strcat(foutput_points, ".points"); + strcat(foutput_faces, ".faces"); + std::ofstream os_points(foutput_points, std::ios::out); + std::ofstream os_faces(foutput_faces, std::ios::out); + if((os_points.fail())||(os_faces.fail())) + std::cerr << "+++unable to open file for output" << std::endl; + else + std::cout << ">> files for output : " << foutput_points + << ", " << foutput_faces << std::endl; + + os_points.clear(); + os_faces.clear(); + + CGAL::set_ascii_mode(os_points); + CGAL::set_ascii_mode(os_faces); + + // af: what is the relationship to _vh_number in Extract_surface + std::size_t _vh_number = std::count_if(T.finite_vertices_begin(), + T.finite_vertices_end(), + Is_not_exterior()); + + os_points << _vh_number << std::endl; + + CGAL::Unique_hash_map vertex_index_map(-1, T.number_of_vertices()); + + int count(0); + for (Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++){ + typename CGAL::Unique_hash_map::Data& d = vertex_index_map[v_it]; + if ((!v_it->is_exterior()) && d == -1){ + d = count; + count++; + os_points << v_it->point() << " 0" << std::endl; + } + } + + os_faces << S.number_of_facets() << std::endl; + for(Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + n = c->neighbor(ci); + ni = n->index(c); + int i1, i2 ,i3; + + if (c->is_selected_facet(ci)) + { + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + os_faces << 3 << " "; + os_faces << vertex_index_map[c->vertex(i1)] + 1 << " "; + os_faces << vertex_index_map[c->vertex(i2)] + 1 << " "; + os_faces << vertex_index_map[c->vertex(i3)] + 1 << " "; + os_faces << " 0 0 0 0" << std::endl; + } + + if (n->is_selected_facet(ni)) + { + i1 = (ni+1) & 3; + i2 = (ni+2) & 3; + i3 = (ni+3) & 3; + os_faces << 3 << " "; + os_faces << vertex_index_map[n->vertex(i1)] + 1 << " "; + os_faces << vertex_index_map[n->vertex(i2)] + 1 << " "; + os_faces << vertex_index_map[n->vertex(i3)] + 1 << " "; + os_faces << " 0 0 0 0" << std::endl; + } + } + + std::cout << "-- medit result written." << std::endl; +} + +//--------------------------------------------------------------------- + +template +void +write_to_file_gv(char* foutput, const Surface& S) +{ + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Vertex Vertex; + typedef typename Surface::Cell_handle Cell_handle; + Triangulation_3& T = S.triangulation(); + + char foutput_tmp[100]; + std::strcpy(foutput_tmp, foutput); + + strcat(foutput_tmp, ".off"); + std::ofstream os(foutput_tmp, std::ios::out); + + if(os.fail()) + std::cerr << "+++unable to open file for output" << std::endl; + else + std::cout << ">> file for output : " << foutput_tmp << std::endl; + + os.clear(); + + CGAL::set_ascii_mode(os); + + std::size_t _vh_number = std::count_if(T.finite_vertices_begin(), + T.finite_vertices_end(), + Is_not_exterior()); + // Header. + os << "OFF" << std::endl + << _vh_number << " " << S.number_of_facets() << " " << 0 << std::endl; + + CGAL::Unique_hash_map vertex_index_map(-1, T.number_of_vertices()); + + int count(0); + for (Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++){ + typename CGAL::Unique_hash_map::Data& d = vertex_index_map[v_it]; + if ((!v_it->is_exterior()) && d == -1){ + d = count; + count++; + os << v_it->point() << " \n"; + } + } + + for(Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + n = c->neighbor(ci); + ni = n->index(c); + int i1, i2 ,i3; + + if (c->is_selected_facet(ci)) + { + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + os << 3 << " "; + os << vertex_index_map[c->vertex(i1)] << " "; + os << vertex_index_map[c->vertex(i2)] << " "; + os << vertex_index_map[c->vertex(i3)] << "\n"; + } + + if (n->is_selected_facet(ni)) + { + i1 = (ni+1) & 3; + i2 = (ni+2) & 3; + i3 = (ni+3) & 3; + os << 3 << " "; + os << vertex_index_map[n->vertex(i1)] << " "; + os << vertex_index_map[n->vertex(i2)] << " "; + os << vertex_index_map[n->vertex(i3)] << "\n"; + } + } + + std::cout << "-- oogl result written." << std::endl; +} + + +template +OutputIterator +write_triple_indices(OutputIterator out, const Surface& S) +{ + std::cerr << "write triple indices" << std::endl; +#if 1 + typedef typename Surface::TDS_2 TDS_2; + typedef typename TDS_2::Face_iterator Face_iterator; + typedef typename Surface::Cell_handle Cell_handle; + + const TDS_2& tds = S.tds_2(); + + for(Face_iterator fit = tds.faces_begin(); fit != tds.faces_end(); ++fit){ + + if(fit->is_on_surface()){ + *out++ = CGAL::Triple(fit->vertex(0)->vertex_3()->id(), + fit->vertex(1)->vertex_3()->id(), + fit->vertex(2)->vertex_3()->id()); + } + } + return out; + +#else + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Vertex Vertex; + typedef typename Surface::Cell_handle Cell_handle; + Triangulation_3& T = S.triangulation(); + + + for(Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + n = c->neighbor(ci); + ni = n->index(c); + int i1, i2 ,i3; + + if (c->is_selected_facet(ci)) + { + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + *out++ = CGAL::Triple(c->vertex(i1)->id(), + c->vertex(i2)->id(), + c->vertex(i3)->id()); + } + + if (n->is_selected_facet(ni)) + { + i1 = (ni+1) & 3; + i2 = (ni+2) & 3; + i3 = (ni+3) & 3; + + *out++ = CGAL::Triple(n->vertex(i1)->id(), + n->vertex(i2)->id(), + n->vertex(i3)->id()); + } + } + return out; +#endif +} + +//--------------------------------------------------------------------- +template +void +write_to_file_ply(char* foutput, const Surface& S) +{ + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Cell_handle Cell_handle; + Triangulation_3& T = S.triangulation(); + char foutput_tmp[100]; + std::strcpy(foutput_tmp, foutput); + + strcat(foutput_tmp, ".ply"); + std::ofstream os(foutput_tmp, std::ios::out | std::ios::binary); + + if(os.fail()) + std::cerr << "+++unable to open file for output" << std::endl; + else + std::cout << ">> file for output : " << foutput_tmp << std::endl; + + os.clear(); + + CGAL::set_ascii_mode(os); + + // Header. + os << "ply" << std::endl + << "format binary_little_endian 1.0" << std::endl + << "comment generated by ply_writer" << std::endl + << "element vertex " << S.number_of_vertices() << std::endl + << "property float x" << std::endl + << "property float y" << std::endl + << "property float z" << std::endl + << "element face " << S.number_of_facets() << std::endl + << "property list uchar int vertex_indices" << std::endl + << "end_header" << std::endl; + + CGAL::set_binary_mode(os); + + CGAL::Unique_hash_map vertex_index_map(-1, T.number_of_vertices()); + + int count(0); + for (Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++){ + typename CGAL::Unique_hash_map::Data& d = vertex_index_map[v_it]; + if ((!v_it->is_exterior()) && d == -1){ + d = count; + count++; + os << v_it->point() << std::endl; + } + } + + for(Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + n = c->neighbor(ci); + ni = n->index(c); + int i1, i2 ,i3; + + if (c->is_selected_facet(ci)) + { + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + char three = '3'; + CGAL::write(os, three, CGAL::io_Read_write()); + CGAL::write(os,vertex_index_map[c->vertex(i1)], CGAL::io_Read_write()); + CGAL::write(os,vertex_index_map[c->vertex(i2)], CGAL::io_Read_write()); + CGAL::write(os,vertex_index_map[c->vertex(i3)], CGAL::io_Read_write()); + os << std::endl; // without color. + // os << 4 << drand48() << drand48() << drand48() << 1.0; // random + // color + } + + if (n->is_selected_facet(ni)) + { + i1 = (ni+1) & 3; + i2 = (ni+2) & 3; + i3 = (ni+3) & 3; + char three = '3'; + CGAL::write(os, three, CGAL::io_Read_write()); + CGAL::write(os,vertex_index_map[n->vertex(i1)], CGAL::io_Read_write()); + CGAL::write(os,vertex_index_map[n->vertex(i2)], CGAL::io_Read_write()); + CGAL::write(os,vertex_index_map[n->vertex(i3)], CGAL::io_Read_write()); + os << std::endl; // without color. + // os << 4 << drand48() << drand48() << drand48() << 1.0; // random + // color + } + } + + //std::cout << "-- ply result written." << std::endl; +} + +//--------------------------------------------------------------------- +template +void +write_to_file_iv_border_edges(const Surface& S, std::ofstream& os) +{ + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_edges_iterator Finite_edges_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Cell_handle Cell_handle; + typedef typename Surface::Edge_like Edge_like; + typedef typename Surface::Border_elt Border_elt; + + Triangulation_3& T = S.triangulation(); + typedef std::pair indiced_vh; + std::map _vh_vect; + int _vh_bord_count = 0; + + for (Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++) + if (v_it->is_on_border()) + { + _vh_vect.insert(indiced_vh (v_it, _vh_bord_count)); + _vh_bord_count++; + } + + typedef const typename Triangulation_3::Point* Const_point_star; + std::vector points_tab(_vh_bord_count); + for (typename std::map::iterator vh_it = _vh_vect.begin(); + vh_it != _vh_vect.end(); vh_it++) + points_tab[vh_it->second] = &vh_it->first->point(); + + os << " Separator {" << std::endl << +" Switch {" << std::endl << +" whichChild 0" << std::endl << +" Separator {" << std::endl << +" BaseColor {" << std::endl << +" rgb 1 0 0" << std::endl << +" }" << std::endl << +" Coordinate3 {" << std::endl << +" point [ "; + bool first(true); + for(int vh_i=0; vh_i<_vh_bord_count; vh_i++) + { + if (!first) os << "," << std::endl << +" "; + else + first=false; + os << *points_tab[vh_i]; + } + os << " ]" << std::endl << +" }" << std::endl << +" IndexedLineSet {" << std::endl << +" coordIndex [ "; + + first=true; + for(Finite_edges_iterator e_it=T.finite_edges_begin(); + e_it!=T.finite_edges_end(); + e_it++) + { + Cell_handle c = (*e_it).first; + int i1 = (*e_it).second, i2 = (*e_it).third; + Edge_like key(c->vertex(i1), c->vertex(i2)); + Border_elt result; + + if (S.is_border_elt(key, result)) + { + if (!first) + os << "," << std::endl << " "; + else + first=false; + os << _vh_vect.find(c->vertex(i1))->second << ", "; + os << _vh_vect.find(c->vertex(i2))->second << ", "; + os << -1; + } + } + os << " ]" << std::endl << +" }}" << std::endl << +" }}" << std::endl; +} + +//--------------------------------------------------------------------- +template +void +write_to_file_iv_remaining_points(const Surface& S, std::ofstream& os) +{ + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Cell_handle Cell_handle; + Triangulation_3& T = S.triangulation(); + typedef std::pair indiced_vh; + std::map _vh_vect; + int _vh_bord_count(0); + + for (Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++) + if (v_it->is_exterior()) + { + _vh_vect.insert(indiced_vh (v_it, _vh_bord_count)); + _vh_bord_count++; + } + + typedef const typename Triangulation_3::Point* Const_point_star; + std::vector points_tab(_vh_bord_count); + for (typename std::map::iterator vh_it = _vh_vect.begin(); + vh_it != _vh_vect.end(); vh_it++) + points_tab[vh_it->second] = &vh_it->first->point(); + + os << " Separator {" << std::endl << +" Switch {" << std::endl << +" whichChild 0" << std::endl << +" Separator {" << std::endl << +" BaseColor {" << std::endl << +" rgb 0 0 1" << std::endl << +" }" << std::endl << +" Coordinate3 {" << std::endl << +" point [ "; + bool first(true); + for(int vh_i=0; vh_i<_vh_bord_count; vh_i++) + { + if (!first) os << "," << std::endl << +" "; + else + first=false; + os << *points_tab[vh_i]; + } + os << " ]" << std::endl << +" }" << std::endl << +" PointSet {" << std::endl << +" startIndex 0" << std::endl << +" numPoints -1" << std::endl << +" }"; + + os << " }" << std::endl << +" }}" << std::endl; + +} + +//--------------------------------------------------------------------- +// attention cette procedure produit un fichier tres sale... trop de sommets... + +// !!!! bizarre : ca a l'air de buggue pour hand.xyz (seg fault...) +template +void +write_to_file_iv_border_facets(const Surface& S, std::ofstream& os) +{ + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Cell_handle Cell_handle; + typedef typename Surface::Edge_like Edge_like; + typedef typename Surface::Border_elt Border_elt; + + Triangulation_3& T = S.triangulation(); + typedef std::pair indiced_vh; + std::map _vh_vect; + int _vh_bord_count(0); + + for (Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++) +// if (v_it->number_of_incident_border() > 0) + { + _vh_vect.insert(indiced_vh (v_it, _vh_bord_count)); + _vh_bord_count++; + } + + typedef const typename Triangulation_3::PoinPoint* Const_point_star; + std::vector points_tab(_vh_bord_count); + for (typename std::map::iterator vh_it = _vh_vect.begin(); + vh_it != _vh_vect.end(); vh_it++) + points_tab[vh_it->second] = &vh_it->first->point(); + + os << " Separator {" << std::endl << +" Switch {" << std::endl << +" whichChild 0" << std::endl << +" Separator {" << std::endl << +" ShapeHints {" << std::endl << +" vertexOrdering CLOCKWISE" << std::endl << +" shapeType UNKNOWN_SHAPE_TYPE" << std::endl << +" faceType CONVEX" << std::endl << +" creaseAngle 1.0" << std::endl << +" }" << std::endl << +" BaseColor {" << std::endl << +" rgb 0 0 1" << std::endl << +" }" << std::endl << +" Coordinate3 {" << std::endl << +" point [ "; + bool first(true); + for(int vh_i=0; vh_i<_vh_bord_count; vh_i++) + { + if (!first) os << "," << std::endl << +" "; + else + first=false; + os << *points_tab[vh_i]; + } + os << " ]" << std::endl << +" }" << std::endl << +" IndexedFaceSet {" << std::endl << +" coordIndex [ "; + + first=true; + for(Finite_facets_iterator f_it=T.finite_facets_begin(); + f_it!=T.finite_facets_end(); + f_it++) + { + Cell_handle c = (*f_it).first; + int index = (*f_it).second; + int i1 = (index+1) & 3; + int i2 = (index+2) & 3; + int i3 = (index+3) & 3; + Edge_like key12(c->vertex(i1), c->vertex(i2)); + Edge_like key13(c->vertex(i1), c->vertex(i3)); + Edge_like key32(c->vertex(i3), c->vertex(i2)); + Border_elt result; + + // les trois aretes sur le bord... +// if (is_border_elt(key12, result)&& +// is_border_elt(key13, result)&& +// is_border_elt(key32, result)) + + // au moins 2 aretes sur le bord... +// if (((is_border_elt(key12, result)&& +// is_border_elt(key13, result)))|| +// ((is_border_elt(key32, result)&& +// is_border_elt(key13, result)))|| +// ((is_border_elt(key12, result)&& +// is_border_elt(key32, result)))) + + // une arete sur le bord... + if ((is_border_elt(key12, result)|| + is_border_elt(key13, result)|| + is_border_elt(key32, result))&& + (c->is_selected_facet(index)|| + c->neighbor(index)->is_selected_facet(c->neighbor(index)->index(c)))) + + // au moins 2 aretes sur le bord... +// if (((is_border_elt(key12, result)&& +// is_border_elt(key13, result)))|| +// ((is_border_elt(key32, result)&& +// is_border_elt(key13, result)))|| +// ((is_border_elt(key12, result)&& +// is_border_elt(key32, result)))) + { + if (!first) + os << "," << std::endl << " "; + else + first=false; + os << _vh_vect.find(c->vertex(i1))->second << ", "; + os << _vh_vect.find(c->vertex(i2))->second << ", "; + os << _vh_vect.find(c->vertex(i3))->second << ", "; + + os << -1; + } + } + os << " ]" << std::endl << +" }}" << std::endl << +" }}" << std::endl; +} + +//--------------------------------------------------------------------- +template +void +write_to_file_iv(char* foutput, const Surface& S, + const bool& boundary) +{ + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Cell_handle Cell_handle; + Triangulation_3& T = S.triangulation(); + char foutput_tmp[100]; + std::strcpy(foutput_tmp, foutput); + + strcat(foutput_tmp, ".iv"); + std::ofstream os(foutput_tmp, std::ios::out); + + if(os.fail()) + std::cerr << "+++unable to open file for output" << std::endl; + else + std::cout << ">> file for output : " << foutput_tmp << std::endl; + + os.clear(); + + CGAL::set_ascii_mode(os); + + // Header. + os << +"#Inventor V2.1 ascii" << std::endl << +"Separator {" << std::endl << +" PerspectiveCamera {" << std::endl << +" position 0 0 2.41421" << std::endl << +" nearDistance 1.41421" << std::endl << +" farDistance 3.41421" << std::endl << +" focalDistance 2.41421" << std::endl << +" }" << std::endl << +" Group {" << std::endl << +" Rotation {" << std::endl << +" }" << std::endl << +" DirectionalLight {" << std::endl << +" direction 0.2 -0.2 -0.979796" << std::endl << +" }" << std::endl << +" ResetTransform {" << std::endl << +" } }" << std::endl << +" Separator {" << std::endl << +" Switch {" << std::endl << +" whichChild 0" << std::endl << +" Separator {" << std::endl << +" ShapeHints {" << std::endl << +" vertexOrdering CLOCKWISE" << std::endl << +" shapeType UNKNOWN_SHAPE_TYPE" << std::endl << +" faceType CONVEX" << std::endl << +" creaseAngle 1.0" << std::endl << +" }" << std::endl << +" BaseColor {" << std::endl << +" rgb 0.6 0.6 0.48" << std::endl << +" }" << std::endl << +" Coordinate3 {" << std::endl << +" point [ "; + + CGAL::Unique_hash_map vertex_index_map(-1, T.number_of_vertices()); + + int count(0); + for (Finite_vertices_iterator v_it = T.finite_vertices_begin(); + v_it != T.finite_vertices_end(); + v_it++){ + typename CGAL::Unique_hash_map::Data& d = vertex_index_map[v_it]; + if ((!v_it->is_exterior()) && d == -1){ + d = count; + count++; + os << v_it->point() << " ,\n"; + } + } + os << " ]" << std::endl << +" }" << std::endl << +" IndexedFaceSet {" << std::endl << +" coordIndex [ "; + + for(Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + n = c->neighbor(ci); + ni = n->index(c); + int i1, i2 ,i3; + + if (c->is_selected_facet(ci)) + { + + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + os << vertex_index_map[c->vertex(i1)] << ", "; + os << vertex_index_map[c->vertex(i2)] << ", "; + os << vertex_index_map[c->vertex(i3)] << ", "; + os << -1; + } + + if (n->is_selected_facet(ni)) + { + i1 = (ni+1) & 3; + i2 = (ni+2) & 3; + i3 = (ni+3) & 3; + os << vertex_index_map[n->vertex(i1)] << ", "; + os << vertex_index_map[n->vertex(i2)] << ", "; + os << vertex_index_map[n->vertex(i3)] << ", "; + os << -1; + } + } + + + os << " ]\n" + " }\n" + " }}\n" + " }\n"; + + if (boundary) + { + // pour visualiser les boundaries restant a la fin... + write_to_file_iv_border_edges(S, os); + + // pour visualiser les facettes eventuellement candidates... + // write_to_file_iv_border_facets(S, os); + + // pour afficher les points non selectionnes, ~bruit??? + // write_to_file_iv_remaining_points(S, os); + } + + os << "}" << std::endl; + + std::cout << "-- Inventor result written." << std::endl; +} + + +template +void +write_boundaries(std::ostream& os, const Surface& S) +{ + + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Boundary_iterator Boundary_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Vertex Vertex; + typedef typename Surface::Cell_handle Cell_handle; + + + CGAL::Random random; + for(Boundary_iterator it = S.boundaries_begin(); + it != S.boundaries_end(); + ++it) { + double blue = random.get_double(0,1); + os << + "Shape {\n" + "appearance Appearance {\n" + "material Material { emissiveColor 1 0 " << blue << "}}\n" + "geometry\n" + "IndexedLineSet {\n" + "coord Coordinate {\n" + "point [ " << std::endl; + unsigned int count = 0; + Vertex_handle first = *it; + do { + os << (*it)->point() << std::endl; + ++it; + count++; + } while(*it != first); + os << "]\n" + "}\n" + "coordIndex [\n"; + + for(unsigned int i = 0; i < count; i++){ + os << i << ", "; + } + os << "0, -1\n"; + os << "]\n" + "}#IndexedLineSet\n" + "}# Shape\n"; + } +} + + + + +template +void +write_to_file_vrml2(char* foutput, const Surface& S, + const bool& boundary, double red, double green, double blue, bool no_header) +{ + + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Cell_handle Cell_handle; + //Triangulation_3& T = S.triangulation(); + + typedef CGAL::Exact_predicates_inexact_constructions_kernel Kernel; + + typedef CGAL::Triangulation_data_structure_2 > TDS; + + TDS tds; + + TDS::Vertex_handle inf = AFSR::orient(tds, S); + + char foutput_tmp[100]; + std::strcpy(foutput_tmp, foutput); + + strcat(foutput_tmp, ".wrl"); + std::ofstream os(foutput_tmp, std::ios::out); + + if(os.fail()) + std::cerr << "+++unable to open file for output" << std::endl; + else + std::cout << ">> file for output : " << foutput_tmp << std::endl; + + os.clear(); + + CGAL::set_ascii_mode(os); + + if(! no_header){ + // Header. + os << + "#VRML V2.0 utf8\n" + "Background {skyColor .1 .5 .5}\n" + "Group {\n" + "children [\n" << std::endl; + }; + + afsr_vrml_output(tds, os, red, green, blue, inf, true); + + if (boundary){ + write_boundaries(os, S); + } + + typename Surface::Outlier_iterator pit = S.outliers_begin(); + + if(pit != S.outliers_end()) { + os << "Shape {\n" + "geometry PointSet {\n" + "coord Coordinate { point [\n"; + + for(; pit != S.outliers_end(); pit++){ + os << pit->x() << " " << pit->y() << " " << pit->z() << ",\n"; + } + os << "] } }\n" + "appearance Appearance {\n" + " material Material {\n" + " emissiveColor 1 0.1 0\n" + " }\n" + "}\n" + "} # Shape\n"; + } + + if(! no_header){ + os << "] # children\n" + "} # Group\n"; + } + std::cout << "-- wrl result written." << std::endl; +} + + + +template +void +write_to_file_stl(char* foutput, const Surface& S) +{ + + typedef typename Surface::Triangulation_3 Triangulation_3; + typedef typename Surface::Finite_vertices_iterator Finite_vertices_iterator; + typedef typename Surface::Finite_facets_iterator Finite_facets_iterator; + typedef typename Surface::Vertex_handle Vertex_handle; + typedef typename Surface::Cell_handle Cell_handle; + typedef typename Triangulation_3::Point Point; + typedef typename CGAL::Kernel_traits::Kernel::Vector_3 Vector; + Triangulation_3& T = S.triangulation(); + + + char foutput_tmp[100]; + std::strcpy(foutput_tmp, foutput); + + strcat(foutput_tmp, ".stl"); + std::ofstream os(foutput_tmp, std::ios::out); + + if(os.fail()) + std::cerr << "+++unable to open file for output" << std::endl; + else + std::cout << ">> file for output : " << foutput_tmp << std::endl; + + os.clear(); + + CGAL::set_ascii_mode(os); + + // Header. + os << "solid" << std::endl; + + for(Finite_facets_iterator f_it = T.finite_facets_begin(); + f_it != T.finite_facets_end(); + f_it++) + { + Cell_handle n, c = (*f_it).first; + int ni, ci = (*f_it).second; + + bool selected = false; + if(c->is_selected_facet(ci)){ + selected = true; + } else { + n = c->neighbor(ci); + ni = n->index(c); + if(n->is_selected_facet(ni)) { + selected = true; + c = n; + ci = ni; + } + } + + if(selected){ + int i1, i2 ,i3; + + i1 = (ci+1) & 3; + i2 = (ci+2) & 3; + i3 = (ci+3) & 3; + + Point p = c->vertex(i1)->point(); + Point q = c->vertex(i2)->point(); + Point r = c->vertex(i3)->point(); + // compute normal + Vector n = CGAL::cross_product( q-p, r-p); + Vector norm = n / sqrt( n * n); + os << "outer loop" << std::endl; + os << "facet normal " << norm << std::endl; + os << "vertex " << p << std::endl; + os << "vertex " << q << std::endl; + os << "vertex " << r << std::endl; + os << "endloop\nendfacet" << std::endl; + } + + } + + os << "endsolid" << std::endl; + + std::cout << "-- stl result written." << std::endl; +} + + +//--------------------------------------------------------------------- +template +void +write_to_file(char* foutput, const Surface& S, + const bool& boundary, const int& out_format, + double red, double green, double blue, bool no_header) +{ + switch(out_format) + { + case -2: + // no output file... + return; + case -1: + write_to_file_iv(foutput, S, boundary); + write_to_file_vrml2(foutput, S, boundary, red, green, blue, no_header); + write_to_file_gv(foutput, S); + write_to_file_medit(foutput, S); + //write_to_file_ply(foutput, S); + return; + case 0: + write_to_file_vrml2(foutput, S, boundary, red, green, blue, no_header); + return; + case 1: + write_to_file_gv(foutput, S); + return; + case 2: + write_to_file_medit(foutput, S); + return; + case 3: + write_to_file_ply(foutput, S); + return; + case 4: + write_to_file_iv(foutput, S, boundary); + return; + case 5: + write_to_file_stl(foutput, S); + return; + } +} + +} // namespace CGAL + + +#endif // CGAL_IO_ADVANCING_FRONT_SURFACE_RECONSTRUCTION_H diff --git a/Advancing_front_surface_reconstruction/include/CGAL/Tvb_3_2.h b/Advancing_front_surface_reconstruction/include/CGAL/Tvb_3_2.h new file mode 100644 index 00000000000..8fbb59d7739 --- /dev/null +++ b/Advancing_front_surface_reconstruction/include/CGAL/Tvb_3_2.h @@ -0,0 +1,94 @@ +// Copyright (c) 1997 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org); you may redistribute it under +// the terms of the Q Public License version 1.0. +// See the file LICENSE.QPL distributed with CGAL. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $Source: /CVSROOT/CGAL/Packages/Triangulation_2/include/CGAL/Triangulation_vertex_base_2.h,v $ +// $Revision$ $Date$ +// $Name: current_submission $ +// +// Author(s) : Mariette Yvinec + + +#ifndef CGAL_TVB_3_2_H +#define CGAL_TVB_3_2_H + +#include +#include + +namespace CGAL { + +template < typename GT, + typename Vb = Triangulation_ds_vertex_base_2<> > +class Tvb_3_2 + : public Vb + +{ + typedef typename Vb::Triangulation_data_structure Tds; +public: + typedef GT Geom_traits; + typedef typename GT::Point_3 Point; + typedef Tds Triangulation_data_structure; + typedef typename Tds::Face_handle Face_handle; + typedef typename Tds::Vertex_handle Vertex_handle; + + template < typename TDS2 > + struct Rebind_TDS { + typedef typename Vb::template Rebind_TDS::Other Vb2; + typedef Tvb_3_2 Other; + }; + +private: + Point _p; + +public: + Tvb_3_2 () : Vb() {} + Tvb_3_2(const Point & p) : Vb(), _p(p) {} + Tvb_3_2(const Point & p, Face_handle f) + : Vb(f), _p(p) {} + Tvb_3_2(Face_handle f) : Vb(f) {} + + void set_point(const Point & p) { _p = p; } + const Point& point() const { return _p; } + + // the non const version of point() is undocument + // but needed to make the point iterator works + // using Lutz projection scheme + Point& point() { return _p; } + + //the following trivial is_valid to allow + // the user of derived face base classes + // to add their own purpose checking + bool is_valid(bool /* verbose */ = false, int /* level */ = 0) const + {return true;} +}; + +template < class GT, class Vb > +std::istream& +operator>>(std::istream &is, Tvb_3_2 &v) + // non combinatorial information. Default = point +{ + return is >> static_cast(v) >> v.point(); +} + +template < class GT, class Vb > +std::ostream& +operator<<(std::ostream &os, const Tvb_3_2 &v) + // non combinatorial information. Default = point +{ + return os << static_cast(v) << v.point(); +} + + + +} //namespace CGAL + +#endif //CGAL_TVB_3_2_H diff --git a/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/copyright b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/copyright new file mode 100644 index 00000000000..434d5723b62 --- /dev/null +++ b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/copyright @@ -0,0 +1 @@ +INRIA Sophia-Antipolis (France) diff --git a/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/description.txt b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/description.txt new file mode 100644 index 00000000000..0d961765618 --- /dev/null +++ b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/description.txt @@ -0,0 +1 @@ +Surface reconstruction using an advancing front method. diff --git a/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/license.txt b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/license.txt new file mode 100644 index 00000000000..8bb8efcb72b --- /dev/null +++ b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/license.txt @@ -0,0 +1 @@ +GPL (v3 or later) diff --git a/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/maintainer b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/maintainer new file mode 100644 index 00000000000..4f45ddcd216 --- /dev/null +++ b/Advancing_front_surface_reconstruction/package_info/Advancing_front_surface_reconstruction/maintainer @@ -0,0 +1 @@ +Andreas Fabri diff --git a/Documentation/doc/Documentation/packages.txt b/Documentation/doc/Documentation/packages.txt index 5dd03dd6775..aa30565b256 100644 --- a/Documentation/doc/Documentation/packages.txt +++ b/Documentation/doc/Documentation/packages.txt @@ -91,6 +91,7 @@ h1 { \package_listing{Surface_reconstruction_points_3} \package_listing{Skin_surface_3} \package_listing{Mesh_3} +\package_listing{Advancing_front_surface_reconstruction} \section PartGeometryProcessing Geometry Processing diff --git a/Polyhedron/demo/Polyhedron/CMakeLists.txt b/Polyhedron/demo/Polyhedron/CMakeLists.txt index 75807ab8bf3..20534d7cd0c 100644 --- a/Polyhedron/demo/Polyhedron/CMakeLists.txt +++ b/Polyhedron/demo/Polyhedron/CMakeLists.txt @@ -28,7 +28,7 @@ foreach(DEP_PKG AABB_tree STL_Extension GraphicsView Surface_mesher Filtered_ker endforeach() # Include this package's headers first -include_directories( BEFORE ./ ./include ../../include ./CGAL_demo) +include_directories( BEFORE ./ ./include ../../include ./CGAL_demo ../../../Advancing_front_surface_reconstruction/include) add_subdirectory( implicit_functions ) @@ -334,6 +334,10 @@ if(CGAL_Qt4_FOUND AND QT4_FOUND AND OPENGL_FOUND AND QGLVIEWER_FOUND) polyhedron_demo_plugin(pca_plugin Polyhedron_demo_pca_plugin) target_link_libraries(pca_plugin scene_polyhedron_item scene_basic_objects) + qt4_wrap_ui( advancing_frontUI_FILES Polyhedron_demo_advancing_front_plugin.ui) + polyhedron_demo_plugin(advancing_front_plugin Polyhedron_demo_advancing_front_plugin Polyhedron_demo_advancing_front_plugin_impl.cpp ${advancing_frontUI_FILES}) + target_link_libraries(advancing_front_plugin scene_polygon_soup_item scene_points_with_normal_item) + if(EIGEN3_FOUND OR TAUCS_FOUND) polyhedron_demo_plugin(parameterization_plugin Polyhedron_demo_parameterization_plugin) target_link_libraries(parameterization_plugin scene_polyhedron_item scene_textured_polyhedron_item ) diff --git a/Polyhedron/demo/Polyhedron/Polyhedron_demo_advancing_front_plugin.cpp b/Polyhedron/demo/Polyhedron/Polyhedron_demo_advancing_front_plugin.cpp new file mode 100755 index 00000000000..d51aee9cd90 --- /dev/null +++ b/Polyhedron/demo/Polyhedron/Polyhedron_demo_advancing_front_plugin.cpp @@ -0,0 +1,121 @@ +#include "config.h" +#include "Scene_points_with_normal_item.h" +#include "Polyhedron_demo_plugin_helper.h" +#include "Polyhedron_demo_plugin_interface.h" +#include + + +#include +#include +#include +#include +#include +#include + +#include "ui_Polyhedron_demo_advancing_front_plugin.h" + +// Poisson reconstruction method: +// Reconstructs a surface mesh from a point set and writes facet indices into polygon soup. +void advancing_front_reconstruct(const Point_set& points, + double sm_perimeter, + double sm_area, + Scene_polygon_soup_item*); + +class Polyhedron_demo_advancing_front_plugin : + public QObject, + public Polyhedron_demo_plugin_helper +{ + Q_OBJECT + Q_INTERFACES(Polyhedron_demo_plugin_interface) + QAction* actionAdvancingFrontReconstruction; + +public: + void init(QMainWindow* mainWindow, Scene_interface* scene_interface) { + + actionAdvancingFrontReconstruction = new QAction(tr("Advancing Front reconstruction"), mainWindow); + actionAdvancingFrontReconstruction->setObjectName("actionAdvancingFrontReconstruction"); + + Polyhedron_demo_plugin_helper::init(mainWindow, scene_interface); + } + + //! Applicate for Point_sets with normals. + bool applicable() const { + return qobject_cast(scene->item(scene->mainSelectionIndex())); + } + + QList actions() const { + return QList() << actionAdvancingFrontReconstruction; + } + +public slots: + void on_actionAdvancingFrontReconstruction_triggered(); +}; // end class Polyhedron_demo_advancing_front_plugin + + +class Polyhedron_demo_advancing_front_plugin_dialog : public QDialog, private Ui::AdvancingFrontDialog +{ + Q_OBJECT + public: + Polyhedron_demo_advancing_front_plugin_dialog(QWidget* /*parent*/ = 0) + { + setupUi(this); + + } + + double trianglePerimeter() const { return m_inputPerimeter->value(); } + double triangleArea() const { return m_inputArea->value(); } +}; + +void Polyhedron_demo_advancing_front_plugin::on_actionAdvancingFrontReconstruction_triggered() +{ + const Scene_interface::Item_id index = scene->mainSelectionIndex(); + + Scene_points_with_normal_item* point_set_item = + qobject_cast(scene->item(index)); + + if(point_set_item) + { + // Gets point set + Point_set* points = point_set_item->point_set(); + if(!points) return; + + // Gets options + Polyhedron_demo_advancing_front_plugin_dialog dialog; + if(!dialog.exec()) + return; + const double sm_perimeter = dialog.trianglePerimeter(); + const double sm_area = dialog.triangleArea(); + + + QApplication::setOverrideCursor(Qt::WaitCursor); + + // Add polyhedron to scene + Scene_polygon_soup_item* new_item = new Scene_polygon_soup_item(); + + for(Point_set::iterator it = points->begin(); it!= points->end(); ++it){ + new_item->new_vertex(it->x(), it->y(), it->z()); + } + + // Reconstruct point set as a polyhedron + advancing_front_reconstruct(*points, sm_perimeter, sm_area, new_item); + + + new_item->setName(tr("%1 Advancing Front (%2 %3)") + .arg(point_set_item->name()) + .arg(sm_perimeter) + .arg(sm_area)); + new_item->setColor(Qt::lightGray); + scene->addItem(new_item); + + // Hide point set + point_set_item->setVisible(false); + scene->itemChanged(index); + + + QApplication::restoreOverrideCursor(); + } +} + +Q_EXPORT_PLUGIN2(Polyhedron_demo_advancing_front_plugin, Polyhedron_demo_advancing_front_plugin) + +#include "Polyhedron_demo_advancing_front_plugin.moc" diff --git a/Polyhedron/demo/Polyhedron/Polyhedron_demo_advancing_front_plugin.ui b/Polyhedron/demo/Polyhedron/Polyhedron_demo_advancing_front_plugin.ui new file mode 100755 index 00000000000..b66834514ab --- /dev/null +++ b/Polyhedron/demo/Polyhedron/Polyhedron_demo_advancing_front_plugin.ui @@ -0,0 +1,113 @@ + + + AdvancingFrontDialog + + + + 0 + 0 + 376 + 170 + + + + Advancing front reconstruction + + + + + + Min triangle perimeter: + + + + + + + * average spacing + + + 0.000000000000000 + + + 30.000000000000000 + + + 0.000000000000000 + + + + + + + Max triangle area: + + + + + + + * average spacing^2 + + + 0 + + + 0.000000000000000 + + + 20.000000000000000 + + + 1.000000000000000 + + + 0.000000000000000 + + + + + + + QDialogButtonBox::Cancel|QDialogButtonBox::Ok + + + + + + + + + buttonBox + accepted() + AdvancingFrontDialog + accept() + + + 177 + 123 + + + 53 + 125 + + + + + buttonBox + rejected() + AdvancingFrontDialog + reject() + + + 257 + 119 + + + 257 + 143 + + + + + diff --git a/Polyhedron/demo/Polyhedron/Scene_points_with_normal_item.cpp b/Polyhedron/demo/Polyhedron/Scene_points_with_normal_item.cpp index 71d86abd0b9..bb7587e934b 100644 --- a/Polyhedron/demo/Polyhedron/Scene_points_with_normal_item.cpp +++ b/Polyhedron/demo/Polyhedron/Scene_points_with_normal_item.cpp @@ -275,11 +275,14 @@ QMenu* Scene_points_with_normal_item::contextMenu() actionResetSelection = menu->addAction(tr("Reset Selection")); actionResetSelection->setObjectName("actionResetSelection"); connect(actionResetSelection, SIGNAL(triggered()),this, SLOT(resetSelection())); +<<<<<<< HEAD actionSelectDuplicatedPoints = menu->addAction(tr("Select duplicated points")); actionSelectDuplicatedPoints->setObjectName("actionSelectDuplicatedPoints"); connect(actionSelectDuplicatedPoints, SIGNAL(triggered()),this, SLOT(selectDuplicates())); +======= +>>>>>>> old menu->setProperty(prop_name, true); } diff --git a/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.cpp b/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.cpp index f86a37b2d54..ccfbf0680a7 100644 --- a/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.cpp +++ b/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.cpp @@ -404,17 +404,19 @@ Scene_polygon_soup_item::bbox() const { } void -Scene_polygon_soup_item::new_vertex(const double& x, - const double& y, - const double& z) +Scene_polygon_soup_item::new_vertex(double x, + double y, + double z) { + if(!soup) + soup = new Polygon_soup; soup->points.push_back(Point_3(x, y, z)); } void -Scene_polygon_soup_item::new_triangle(const std::size_t i, - const std::size_t j, - const std::size_t k) +Scene_polygon_soup_item::new_triangle(std::size_t i, + std::size_t j, + std::size_t k) { Polygon_soup::Polygon_3 new_polygon(3); new_polygon[0] = i; diff --git a/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.h b/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.h index 648056d1e3f..62fe4139ff7 100644 --- a/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.h +++ b/Polyhedron/demo/Polyhedron/Scene_polygon_soup_item.h @@ -34,10 +34,10 @@ public: bool isFinite() const { return true; } bool isEmpty() const; Bbox bbox() const; - - void new_vertex(const double&, const double&, const double&); - void new_triangle(const std::size_t, const std::size_t, const std::size_t); - + + void new_vertex(double, double, double); + void new_triangle(std::size_t, std::size_t, std::size_t); + public slots: void shuffle_orientations(); bool orient(); diff --git a/Triangulation_2/include/CGAL/Triangulation_data_structure_2.h b/Triangulation_2/include/CGAL/Triangulation_data_structure_2.h index 63ea3e99578..a0a9f788e86 100644 --- a/Triangulation_2/include/CGAL/Triangulation_data_structure_2.h +++ b/Triangulation_2/include/CGAL/Triangulation_data_structure_2.h @@ -301,7 +301,9 @@ public: Face_handle f1, Face_handle f2, Face_handle f3); + void set_adjacency(Face_handle f0, int i0, Face_handle f1, int i1) const; + void delete_face(Face_handle); void delete_vertex(Vertex_handle); @@ -342,10 +344,12 @@ public: // HELPING private: typedef std::pair Vh_pair; +public: void set_adjacency(Face_handle fh, int ih, std::map< Vh_pair, Edge>& edge_map); void reorient_faces(); +private: bool dim_down_precondition(Face_handle f, int i); public: @@ -2224,24 +2228,24 @@ reorient_faces() std::set oriented_set; std::stack st; Face_iterator fit = faces_begin(); - int nf = std::distance(faces_begin(),faces_end()); + std::ptrdiff_t nf = std::distance(faces_begin(),faces_end()); - while (static_cast(oriented_set.size()) != nf) { - while ( oriented_set.find(fit) != oriented_set.end()){ + while (0 != nf) { + while ( !oriented_set.insert(fit).second ){ ++fit; // find a germ for non oriented components } // orient component - oriented_set.insert(fit); + --nf; st.push(fit); while ( ! st.empty()) { Face_handle fh = st.top(); st.pop(); for(int ih = 0 ; ih < 3 ; ++ih){ Face_handle fn = fh->neighbor(ih); - if (oriented_set.find(fn) == oriented_set.end()){ + if (oriented_set.insert(fn).second){ int in = fn->index(fh); if (fn->vertex(cw(in)) != fh->vertex(ccw(ih))) fn->reorient(); - oriented_set.insert(fn); + --nf; st.push(fn); } }