40#ifndef GEOGRAM_MESH_MESH_SURFACE_INTERSECTION
41#define GEOGRAM_MESH_MESH_SURFACE_INTERSECTION
49#include <geogram/basic/debug_stream.h>
191 monster_threshold_ = nb;
219 detect_intersecting_neighbors_ = x;
230 use_radial_sort_ = x;
255 skeleton_ = skeleton;
256 skeleton_trim_fins_ = trim_fins;
265 interpolate_attributes_ = x;
340 Process::acquire_spinlock(lock_);
347 Process::release_spinlock(lock_);
513 Sign o_ref1 = h_orient(h_ref_, h1);
514 Sign o_ref2 = h_orient(h_ref_, h2);
515 Sign oN_ref1 = h_refNorient(h1);
516 Sign oN_ref2 = h_refNorient(h2);
517 Sign o_12 = h_orient(h1,h2);
519 <<
" o_ref1=" << int(o_ref1) <<
" o_ref2=" << int(o_ref2)
520 <<
" oN_ref1=" << int(oN_ref1) <<
" oN_ref2=" << int(oN_ref2)
521 <<
" o_12=" << int(o_12)
536 mutable bool degenerate_;
546#ifdef GEOGRAM_USE_EXACT_NT
551 typedef vec3HExLexicoCompareCanonical ExactPointCompare;
556 std::map<ExactPoint,index_t,ExactPointCompare> exact_point_to_vertex_;
561 bool detect_intersecting_neighbors_;
562 bool use_radial_sort_;
567 vec3 normalize_center_;
568 double normalize_radius_;
576 bool skeleton_trim_fins_;
577 bool interpolate_attributes_;
611 facet_corner_alpha3_.bind(
622 return mesh_.facet_corners.
nb();
683 return facet_corner_alpha3_[h];
693 return alpha3(3*f)/3;
709 return mesh_.facets.
vertex(f,lv);
741 facet_corner_alpha3_[h1] = h2;
742 facet_corner_alpha3_[h2] = h1;
777 return bndl_start_.size() - 1;
803 return bndl_start_[bndl+1] - bndl_start_[bndl];
816 return H_[bndl_start_[bndl] + li];
829 H_[bndl_start_[bndl] + li] = h;
839 H_, bndl_start_[bndl], bndl_start_[bndl+1]
850 H_, bndl_start_[bndl], bndl_start_[bndl+1]
863 index_t h = H_[bndl_start_[bndl]];
864 return I_.halfedges_.vertex(h,lv);
875 return v_first_bndl_[v];
886 return bndl_next_around_v_[bndl];
897 index_t bndl = vertex_first_bundle(v);
899 bndl = next_around_vertex(bndl)
914 return (bndl >= nb()/2) ? (bndl-nb()/2) : (bndl+nb()/2);
924 if(nb_bundles_around_vertex(v) != 2) {
928 index_t bndl2 = vertex_first_bundle(v);
929 bndl2 != NO_INDEX; bndl2 = next_around_vertex(bndl2)
932 return opposite(bndl2);
945 if(nb_bundles_around_vertex(v) != 2) {
949 index_t bndl2 = vertex_first_bundle(v);
950 bndl2 != NO_INDEX; bndl2 = next_around_vertex(bndl2)
952 if(opposite(bndl2) != bndl) {
969 if(nb_halfedges(bndl) <= 2) {
972 auto b = H_.begin() + std::ptrdiff_t(bndl_start_[bndl]);
973 auto e = H_.begin() + std::ptrdiff_t(bndl_start_[bndl+1]);
982 bndl_is_sorted_[bndl] = OK;
999 set_halfedge(bndl, i, halfedges[i]);
1001 bndl_is_sorted_[bndl] =
true;
1022 bool is_sorted(
index_t bndl)
const {
1024 return bndl_is_sorted_[bndl];
1066 return polyline_start_.size() - 1;
1093 B_, polyline_start_[polyline], polyline_start_[polyline+1]
1099 return polyline_start_[polyline+1] - polyline_start_[polyline];
1105 return B_[polyline_start_[polyline] + li];
1122 } radial_polylines_;
1127 enum MeshBooleanOperationFlags {
1128 MESH_BOOL_OPS_DEFAULT = 0,
1129 MESH_BOOL_OPS_VERBOSE = 1,
1130 MESH_BOOL_OPS_ATTRIBS = 2,
1131 MESH_BOOL_OPS_NO_SIMPLIFY = 4
1148 Mesh& result,
Mesh& A,
Mesh& B,
const std::string& operation,
1149 MeshBooleanOperationFlags flags = MESH_BOOL_OPS_DEFAULT
1163 Mesh& result,
Mesh& A,
Mesh& B,
const std::string& operation,
1167 result, A, B, operation,
1168 verbose ? MESH_BOOL_OPS_VERBOSE : MESH_BOOL_OPS_DEFAULT
1186 MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT
1221 MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT
1255 MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT
#define geo_assert_not_reached
Sets a non reachable point in the program.
#define geo_assert(x)
Verifies that a condition is met.
#define geo_debug_assert(x)
Verifies that a condition is met.
Generic mechanism for attributes.
Common include file, providing basic definitions. Should be included before anything else by all head...
Manages an attribute attached to a set of object.
Detects and retriangulates a set of coplanar facets for MeshSurfaceIntersection.
index_t vertex(index_t f, index_t lv) const
Gets a vertex by facet and local vertex index.
index_range corners(index_t f) const
Gets the corners of a facet.
Meshes a single triangle with the constraints that come from the intersections with the other triangl...
AttributesManager & attributes() const
Gets the attributes manager.
index_t nb() const
Gets the number of (sub-)elements.
Halfedfge-like API wrappers on top of a triangulated mesh.
index_t facet(index_t h) const
Gets the facet associated to a halfedge.
void sew2(index_t h1, index_t h2)
Creates a surfacic link between two halfedges.
void sew3(index_t h1, index_t h2)
Creates a volumetric link between two halfedges.
void initialize()
Initializes the structure.
index_t alpha3(index_t h) const
gets the volumetric neighbor of a halfedge
index_t vertex(index_t h, index_t dlv) const
gets a vertex of an halfedge
index_t alpha2(index_t h) const
gets the surfacic neighbor of a halfedge
index_t facet_alpha3(index_t f) const
gets the volumetric neighbor of a facet
Halfedges(MeshSurfaceIntersection &I)
Halfedges constructor.
index_as_iterator begin() const
used by range-based for
index_t nb() const
Gets the number of halfedegs in the map.
~Halfedges()
Halfedges destructor.
index_as_iterator end() const
used by range-based for
Represents the set of radial halfedge bundles.
index_as_iterator begin() const
used by range-based for
index_t nb() const
Gets the number of bundles.
index_t vertex(index_t bndl, index_t lv) const
gets one of the vertices at the two extremities of a bundle
index_t vertex_first_bundle(index_t v) const
gets the first bundle starting from a vertex
void get_sorted_incident_charts(index_t bndl, vector< ChartPos > &chart_pos)
Gets the sorted list of charts around bundle.
index_as_iterator end() const
used by range-based for
const_index_ptr_range halfedges(index_t bndl) const
gets the halfedges in a bundle
index_t next_around_vertex(index_t bndl) const
gets the next bundle around a vertex
bool radial_sort(index_t bndl, RadialSort &RS)
Sorts the halfedges of the bundle in-place.
index_ptr_range halfedges(index_t bndl)
gets the halfedges in a bundle
index_t prev_along_polyline(index_t bndl)
gets the predecessor of a bundle along its polyline
index_t nb_halfedges(index_t bndl) const
Gets the number of halfedges in a bundle.
void initialize()
Initializes the structure.
index_t next_along_polyline(index_t bndl)
gets the successor of a bundle along its polyline
void set_sorted_halfedges(index_t bndl, const vector< index_t > &halfedges)
Sets the halfedges of a bundle.
std::pair< index_t, index_t > ChartPos
Indicates where to find a chart in a bundle.
index_t nb_bundles_around_vertex(index_t v) const
gets the bumber of bundles around a vertex
RadialBundles(MeshSurfaceIntersection &I)
RadialBundles constructor.
index_t opposite(index_t bndl)
gets the opposite bundle
index_t halfedge(index_t bndl, index_t li) const
Gets a halfedge in a bundle from local index.
void set_halfedge(index_t bndl, index_t li, index_t h)
Sets a halfedge in a bundle.
void radial_sort()
Sorts all the bundles of all polylines.
index_as_iterator begin() const
used by range-based for
const_index_ptr_range bundles(index_t polyline) const
gets the bundles in a polyline
void get_skeleton(Mesh &to, bool trim_fins=false)
Copies the set of polylines to a mesh.
void initialize()
Initializes the structure.
index_t nb() const
Gets the number of polylines.
index_as_iterator end() const
used by range-based for
RadialPolylines(MeshSurfaceIntersection &I)
RadialPolylines constructor.
Sign h_orient(index_t h1, index_t h2) const
Computes the relative orientations of two halfedges.
static exact::vec3 exact_direction(const ExactPoint &p1, const ExactPoint &p2)
Computes a vector of arbitrary length with its direction given by two points.
Sign h_refNorient(index_t h2) const
Computes the normal orientation of a halfedge relative to h_ref.
void init(index_t h_ref)
Initializes radial sorting around a given halfedge.
bool degenerate() const
Tests if a degeneracy was encountered.
RadialSort(const MeshSurfaceIntersection &mesh)
RadialSort constructor.
bool operator()(index_t h1, index_t h2) const
Compares two halfedges.
static vec3I exact_direction_I(const ExactPoint &p1, const ExactPoint &p2)
Computes an interval vector of arbitrary length with its direction given by two points.
Computes surface intersections.
void set_verbose(bool x)
Display information while computing the intersection. Default is unset.
index_t tentatively_classify_component_vertex_fast(index_t component, index_t v)
Classifies a vertex of the computed intersection.
void set_monster_threshold(index_t nb)
Sets the threshold from which triangle is considered to be a monster.
void simplify_coplanar_facets(double angle_tolerance=0.0)
Merge coplanar facets and retriangulate them using a Constrained Delaunay triangulation.
void remove_internal_shells()
Removes all the facets that are not on the outer boundary.
index_t classify_component(index_t component, index_t v)
Classifies a connected component.
void unlock()
Releases the lock associated with this mesh.
void set_detect_intersecting_neighbors(bool x)
detect and compute intersections between facets that share a facet or an edge. Set to false if input ...
void set_dry_run(bool x)
In dry run mode, the computed local triangulations are not inserted in the global mesh....
index_t find_or_create_exact_vertex(const ExactPoint &p)
Finds or creates a vertex in the mesh, by exact coordinates.
void set_normalize(bool x)
Specifies whether coordinates should be normalized during computation. If set, original coordinates a...
void intersect_prologue()
substep of intersect(), prepares the mesh
void remove_external_shell()
Removes all the facets that are on the outer boundary.
const Mesh & readonly_mesh() const
Gets a copy of the initial mesh passed to the constructor.
void set_build_skeleton(Mesh *skeleton, bool trim_fins=false)
Optionally save the skeleton (that is, the collection of non-manifold edges) to a given mesh....
const Mesh & target_mesh() const
Gets the target mesh.
ExactPoint exact_vertex(index_t v) const
Gets the exact point associated with a vertex.
void set_radial_sort(bool x)
Specifies whether surfaces should be duplicated and radial edges sorted in order to create the volume...
void intersect_epilogue(const vector< IsectInfo > &intersections)
subset of intersect(), cleans the resulting mesh and undoes optional geometric normalization.
void intersect_get_intersections(vector< IsectInfo > &intersections)
substep of intersect(), finds all the intersection points and segments.
void classify(const std::string &expr)
Classifies the facets and keep only the ones on the boundary of a combination of regions defined by a...
index_t tentatively_classify_component_vertex(index_t component, index_t v)
Classifies a vertex of the computed intersection.
void lock()
Acquires a lock on this mesh.
void intersect_remesh_intersections(vector< IsectInfo > &intersections)
substep of intersect(), inserts the intersection points and segments into the triangles.
void mark_external_shell(vector< index_t > &on_external_shell)
Marks all the facets that are on the external shell.
void build_Weiler_model()
Builds the Weiler model.
void set_interpolate_attributes(bool x)
Specifies that attributes should be interpolated.
Mesh & target_mesh()
Gets the target mesh.
void set_delaunay(bool x)
If set, compute constrained Delaunay triangulation in the intersected triangles. If there are interse...
Wraps an integer to be used with the range-based for construct.
Comparator class for vec3Hg \detail Used to create maps indexed by vec3Hg or SOS symbolic perturbatio...
3d vector with homogeneous coordinates
Vector with aligned memory allocation.
index_t size() const
Gets the number of elements.
Exact predicates and constructs.
The class that represents a mesh.
Functions to load and save meshes.
SOSMode
Mode for symbolic perturbations.
std::atomic_flag spinlock
A lightweight synchronization structure.
Global Vorpaline namespace.
bool mesh_facets_have_intersection(Mesh &M, index_t f1, index_t f2)
Tests whether two mesh facets have a non-degenerate intersection.
void mesh_boolean_operation(Mesh &result, Mesh &A, Mesh &B, const std::string &operation, MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT)
Computes a boolean operation with two surface meshes.
void mesh_union(Mesh &result, Mesh &A, Mesh &B, MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT)
Computes the union of two surface meshes.
geo_index_t index_t
The type for storing and manipulating indices.
Sign
Integer constants that represent the sign of a value.
void mesh_remove_intersections(Mesh &M, index_t max_iter=3, bool verbose=false)
Attempts to make a surface mesh conformal by removing intersecting facets and re-triangulating the ho...
void mesh_difference(Mesh &result, Mesh &A, Mesh &B, MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT)
Computes the difference of two surface meshes.
void mesh_intersection(Mesh &result, Mesh &A, Mesh &B, MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT)
Computes the intersection of two surface meshes.
Function and classes for process manipulation.