Geogram Version 1.9.6
A programming library of geometric algorithms
Loading...
Searching...
No Matches
mesh_surface_intersection.h
Go to the documentation of this file.
1/*
2 * Copyright (c) 2000-2022 Inria
3 * All rights reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 *
8 * * Redistributions of source code must retain the above copyright notice,
9 * this list of conditions and the following disclaimer.
10 * * Redistributions in binary form must reproduce the above copyright notice,
11 * this list of conditions and the following disclaimer in the documentation
12 * and/or other materials provided with the distribution.
13 * * Neither the name of the ALICE Project-Team nor the names of its
14 * contributors may be used to endorse or promote products derived from this
15 * software without specific prior written permission.
16 *
17 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
18 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
19 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
20 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
21 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
22 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
23 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
24 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
25 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
26 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
27 * POSSIBILITY OF SUCH DAMAGE.
28 *
29 * Contact: Bruno Levy
30 *
31 * https://www.inria.fr/fr/bruno-levy
32 *
33 * Inria,
34 * Domaine de Voluceau,
35 * 78150 Le Chesnay - Rocquencourt
36 * FRANCE
37 *
38 */
39
40#ifndef GEOGRAM_MESH_MESH_SURFACE_INTERSECTION
41#define GEOGRAM_MESH_MESH_SURFACE_INTERSECTION
42
44#include <geogram/mesh/mesh.h>
49#include <geogram/basic/debug_stream.h>
50#include <functional>
51
58namespace GEO {
59
60 struct IsectInfo;
61
62 /********************************************************************/
63
68 class GEOGRAM_API MeshSurfaceIntersection {
69 public:
70
72
75
83 void intersect();
84
90
96
97
98 void remove_fins();
99
123 void classify(const std::string& expr);
124
125
134
148 index_t component, index_t v
149 );
150
162 index_t component, index_t v
163 );
164
165
172 void simplify_coplanar_facets(double angle_tolerance = 0.0);
173
178 void set_verbose(bool x) {
179 verbose_ = x;
180 if(!verbose_ && fine_verbose_) {
181 fine_verbose_ = false;
182 }
183 }
184
189 void set_fine_verbose(bool x) {
190 fine_verbose_ = x;
191 if(fine_verbose_ && !verbose_) {
192 verbose_ = true;
193 }
194 }
195
204 monster_threshold_ = nb;
205 }
206
212 void set_dry_run(bool x) {
213 dry_run_ = x;
214 }
215
222 void set_delaunay(bool x) {
223 delaunay_ = x;
224 }
225
232 detect_intersecting_neighbors_ = x;
233 }
234
242 void set_radial_sort(bool x) {
243 use_radial_sort_ = x;
244 }
245
253 void set_normalize(bool x) {
254 normalize_ = x;
255 }
256
257
267 void set_build_skeleton(Mesh* skeleton, bool trim_fins=false) {
268 skeleton_ = skeleton;
269 skeleton_trim_fins_ = trim_fins;
270 }
271
278 interpolate_attributes_ = x;
279 }
280
281 protected:
291
306
325
341 void intersect_epilogue(const vector<IsectInfo>& intersections);
342
343
352 void lock() {
353 Process::acquire_spinlock(lock_);
354 }
355
359 void unlock() {
360 Process::release_spinlock(lock_);
361 }
362
373
385
392 return mesh_;
393 }
394
400 const Mesh& target_mesh() const {
401 return mesh_;
402 }
403
416 const Mesh& readonly_mesh() const {
417 return mesh_copy_;
418 }
419
420 class RadialSort;
421
432
436 void mark_external_shell(vector<index_t>& on_external_shell);
437
438 protected:
439
443 class GEOGRAM_API RadialSort {
444 public:
450 const MeshSurfaceIntersection& mesh
451 ) : mesh_(mesh),
452 h_ref_(NO_INDEX),
453 degenerate_(false)
454 {
455 }
456
461 void init(index_t h_ref);
462
469 bool operator()(index_t h1, index_t h2) const;
470
477 bool degenerate() const {
478 return degenerate_;
479 }
480
489 const ExactPoint& p1, const ExactPoint& p2
490 );
491
500 const ExactPoint& p1, const ExactPoint& p2
501 );
502
503 protected:
504
514
522
523 public:
524 void test(index_t h1, index_t h2) {
525 (*this)(h1,h2);
526 Sign o_ref1 = h_orient(h_ref_, h1);
527 Sign o_ref2 = h_orient(h_ref_, h2);
528 Sign oN_ref1 = h_refNorient(h1);
529 Sign oN_ref2 = h_refNorient(h2);
530 Sign o_12 = h_orient(h1,h2);
531 std::cerr
532 << " o_ref1=" << int(o_ref1) << " o_ref2=" << int(o_ref2)
533 << " oN_ref1=" << int(oN_ref1) << " oN_ref2=" << int(oN_ref2)
534 << " o_12=" << int(o_12)
535 << std::endl;
536 }
537
538
539 private:
540 const MeshSurfaceIntersection& mesh_;
541 index_t h_ref_; // ---reference halfedge
542 exact::vec3 U_ref_; // -.
543 exact::vec3 V_ref_; // +-reference basis
544 exact::vec3 N_ref_; // -'
545 vec3I U_ref_I_; // -.
546 vec3I V_ref_I_; // +-reference basis (interval arithmetics)
547 vec3I N_ref_I_; // _'
548 mutable vector< std::pair<index_t, Sign> > refNorient_cache_;
549 mutable bool degenerate_;
550 };
551
552
553 protected:
554 Process::spinlock lock_;
555 Mesh& mesh_;
556 Mesh mesh_copy_;
557 Attribute<const ExactPoint*> vertex_to_exact_point_;
558
559#ifdef GEOGRAM_USE_EXACT_NT
560 // Exact points are canonicalized
561 // (by Numeric::optimize_number_representation(vec3HEx)) so
562 // we can use this comparator that makes the global vertex map
563 // much much faster.
564 typedef vec3HExLexicoCompareCanonical ExactPointCompare;
565#else
566 // Generic comparator for global vertex map.
567 typedef vec3HgLexicoCompare<exact::scalar> ExactPointCompare;
568#endif
569 std::map<ExactPoint,index_t,ExactPointCompare> exact_point_to_vertex_;
570
571 bool verbose_;
572 bool fine_verbose_;
573 bool delaunay_;
574 bool detect_intersecting_neighbors_;
575 bool use_radial_sort_;
576
577 PCK::SOSMode SOS_bkp_;
578 bool rescale_;
579 bool normalize_;
580 vec3 normalize_center_;
581 double normalize_radius_;
582
583 index_t monster_threshold_;
584 bool dry_run_;
585 friend class MeshInTriangle;
586 friend class CoplanarFacets;
587
588 Mesh* skeleton_;
589 bool skeleton_trim_fins_;
590 bool interpolate_attributes_;
591
592 /***************************************************/
593
602 class Halfedges {
603 public:
604
609 Halfedges(MeshSurfaceIntersection& I) : mesh_(I.mesh_) {
610 }
611
616 // TODO: destroy alpha3 attribute (kept now for debugging
617 }
618
623 void initialize() {
624 facet_corner_alpha3_.bind(
625 mesh_.facet_corners.attributes(), "alpha3"
626 );
627 }
628
634 index_t nb() const {
635 return mesh_.facet_corners.nb();
636 }
637
643 return index_as_iterator(0);
644 }
645
651 return index_as_iterator(nb());
652 }
653
660 return h/3;
661 }
662
673 index_t t1 = facet(h);
674 index_t t2 = mesh_.facet_corners.adjacent_facet(h);
675 if(t2 == NO_INDEX) {
676 return NO_INDEX;
677 }
678 for(index_t h2: mesh_.facets.corners(t2)) {
679 if(mesh_.facet_corners.adjacent_facet(h2) == t1) {
680 return h2;
681 }
682 }
684 }
685
696 return facet_corner_alpha3_[h];
697 }
698
706 return alpha3(3*f)/3;
707 }
708
720 index_t f = h/3;
721 index_t lv = (h+dlv)%3;
722 return mesh_.facets.vertex(f,lv);
723 }
724
725
734 void sew2(index_t h1, index_t h2) {
735 geo_debug_assert(vertex(h1,0) == vertex(h2,1));
736 geo_debug_assert(vertex(h2,0) == vertex(h1,1));
737 index_t t1 = h1/3;
738 index_t t2 = h2/3;
739 mesh_.facet_corners.set_adjacent_facet(h1,t2);
740 mesh_.facet_corners.set_adjacent_facet(h2,t1);
741 }
742
751 void sew3(index_t h1, index_t h2) {
752 geo_debug_assert(vertex(h1,0) == vertex(h2,1));
753 geo_debug_assert(vertex(h2,0) == vertex(h1,1));
754 facet_corner_alpha3_[h1] = h2;
755 facet_corner_alpha3_[h2] = h1;
756 }
757
758 private:
759 Mesh& mesh_;
760 Attribute<index_t> facet_corner_alpha3_;
761 } halfedges_;
762
763 /***************************************************/
764
771 public:
772
777 RadialBundles(MeshSurfaceIntersection& I) : I_(I), mesh_(I.mesh_) {
778 }
779
785
789 index_t nb() const {
790 return bndl_start_.size() - 1;
791 }
792
798 return index_as_iterator(0);
799 }
800
806 return index_as_iterator(nb());
807 }
808
815 geo_debug_assert(bndl < nb());
816 return bndl_start_[bndl+1] - bndl_start_[bndl];
817 }
818
826 index_t halfedge(index_t bndl, index_t li) const {
827 geo_debug_assert(bndl < nb());
828 geo_debug_assert(li < nb_halfedges(bndl));
829 return H_[bndl_start_[bndl] + li];
830 }
831
840 geo_debug_assert(bndl < nb());
841 geo_debug_assert(li < nb_halfedges(bndl));
842 H_[bndl_start_[bndl] + li] = h;
843 }
844
851 return index_ptr_range(
852 H_, bndl_start_[bndl], bndl_start_[bndl+1]
853 );
854 }
855
863 H_, bndl_start_[bndl], bndl_start_[bndl+1]
864 );
865 }
866
874 index_t vertex(index_t bndl, index_t lv) const {
875 geo_debug_assert(bndl_start_[bndl+1] - bndl_start_[bndl] > 0);
876 index_t h = H_[bndl_start_[bndl]];
877 return I_.halfedges_.vertex(h,lv);
878 }
879
888 return v_first_bndl_[v];
889 }
890
899 return bndl_next_around_v_[bndl];
900 }
901
908 index_t result = 0;
909 for(
910 index_t bndl = vertex_first_bundle(v);
911 bndl != NO_INDEX;
912 bndl = next_around_vertex(bndl)
913 ) {
914 ++result;
915 }
916 return result;
917 }
918
926 geo_debug_assert(bndl < nb());
927 return (bndl >= nb()/2) ? (bndl-nb()/2) : (bndl+nb()/2);
928 }
929
936 index_t v = vertex(bndl,0);
937 if(nb_bundles_around_vertex(v) != 2) {
938 return NO_INDEX;
939 }
940 for(
941 index_t bndl2 = vertex_first_bundle(v);
942 bndl2 != NO_INDEX; bndl2 = next_around_vertex(bndl2)
943 ) {
944 if(bndl2 != bndl) {
945 return opposite(bndl2);
946 }
947 }
949 }
950
957 index_t v = vertex(bndl,1);
958 if(nb_bundles_around_vertex(v) != 2) {
959 return NO_INDEX;
960 }
961 for(
962 index_t bndl2 = vertex_first_bundle(v);
963 bndl2 != NO_INDEX; bndl2 = next_around_vertex(bndl2)
964 ) {
965 if(opposite(bndl2) != bndl) {
966 return bndl2;
967 }
968 }
970 }
971
981 geo_debug_assert(bndl < nb());
982 if(nb_halfedges(bndl) <= 2) {
983 return true;
984 }
985 auto b = H_.begin() + std::ptrdiff_t(bndl_start_[bndl]);
986 auto e = H_.begin() + std::ptrdiff_t(bndl_start_[bndl+1]);
987 RS.init(*b);
988 std::sort(
989 b, e,
990 [&](index_t h1, index_t h2)->bool {
991 return RS(h1,h2);
992 }
993 );
994 bool OK = !RS.degenerate();
995 bndl_is_sorted_[bndl] = OK;
996 return OK;
997 }
998
1008 index_t bndl, const vector<index_t>& halfedges
1009 ) {
1010 geo_debug_assert(halfedges.size() == nb_halfedges(bndl));
1011 for(index_t i=0; i<halfedges.size(); ++i) {
1012 set_halfedge(bndl, i, halfedges[i]);
1013 }
1014 bndl_is_sorted_[bndl] = true;
1015 }
1016
1022 typedef std::pair<index_t, index_t> ChartPos;
1023
1032 index_t bndl, vector<ChartPos>& chart_pos
1033 );
1034
1035 bool is_sorted(index_t bndl) const {
1036 geo_assert(bndl < nb());
1037 return bndl_is_sorted_[bndl];
1038 }
1039
1040 // private:
1042 Mesh& mesh_;
1043 Attribute<index_t> facet_chart_;
1044 vector<index_t> H_;
1045 vector<index_t> bndl_start_;
1046 vector<index_t> v_first_bndl_;
1047 vector<index_t> bndl_next_around_v_;
1048 vector<bool> bndl_is_sorted_;
1049 } radial_bundles_;
1050
1051 /***************************************************/
1052
1054 public:
1059 RadialPolylines(MeshSurfaceIntersection& I) : I_(I), mesh_(I.mesh_) {
1060 }
1061
1067
1074
1078 index_t nb() const {
1079 return polyline_start_.size() - 1;
1080 }
1081
1087 return index_as_iterator(0);
1088 }
1089
1095 return index_as_iterator(nb());
1096 }
1097
1104 geo_debug_assert(polyline < nb());
1105 return const_index_ptr_range(
1106 B_, polyline_start_[polyline], polyline_start_[polyline+1]
1107 );
1108 }
1109
1110 index_t nb_bundles(index_t polyline) const {
1111 geo_debug_assert(polyline < nb());
1112 return polyline_start_[polyline+1] - polyline_start_[polyline];
1113 }
1114
1115 index_t bundle(index_t polyline, index_t li) const {
1116 geo_debug_assert(polyline < nb());
1117 geo_debug_assert(li < nb_bundles(polyline));
1118 return B_[polyline_start_[polyline] + li];
1119 }
1120
1128 void get_skeleton(Mesh& to, bool trim_fins=false);
1129
1130 private:
1132 Mesh& mesh_;
1133 vector<index_t> B_;
1134 vector<index_t> polyline_start_;
1135 } radial_polylines_;
1136 };
1137
1138 /********************************************************************/
1139
1140 enum MeshBooleanOperationFlags {
1141 MESH_BOOL_OPS_DEFAULT = 0,
1142 MESH_BOOL_OPS_VERBOSE = 1,
1143 MESH_BOOL_OPS_ATTRIBS = 2,
1144 MESH_BOOL_OPS_NO_SIMPLIFY = 4
1145 };
1146
1160 void GEOGRAM_API mesh_boolean_operation(
1161 Mesh& result, const Mesh& A, const Mesh& B, const std::string& operation,
1162 MeshBooleanOperationFlags flags = MESH_BOOL_OPS_DEFAULT
1163 );
1164
1176 Mesh& result, const Mesh& A, const Mesh& B, const std::string& operation,
1177 bool verbose
1178 ) {
1180 result, A, B, operation,
1181 verbose ? MESH_BOOL_OPS_VERBOSE : MESH_BOOL_OPS_DEFAULT
1182 );
1183 }
1184
1197 inline void mesh_union(
1198 Mesh& result, const Mesh& A, const Mesh& B,
1199 MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT
1200 ) {
1201 mesh_boolean_operation(result, A, B, "A+B", flags);
1202 }
1203
1213 inline void mesh_union(
1214 Mesh& result, const Mesh& A, const Mesh& B, bool verbose
1215 ) {
1216 mesh_boolean_operation(result, A, B, "A+B", verbose);
1217 }
1218
1219
1233 Mesh& result, const Mesh& A, const Mesh& B,
1234 MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT
1235 ) {
1236 mesh_boolean_operation(result, A, B, "A*B", flags);
1237 }
1238
1249 Mesh& result, const Mesh& A, const Mesh& B, bool verbose
1250 ) {
1251 mesh_boolean_operation(result, A, B, "A*B", verbose);
1252 }
1253
1266 inline void mesh_difference(
1267 Mesh& result, const Mesh& A, const Mesh& B,
1268 MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT
1269 ) {
1270 mesh_boolean_operation(result, A, B, "A-B", flags);
1271 }
1272
1282 inline void mesh_difference(
1283 Mesh& result, const Mesh& A, const Mesh& B, bool verbose
1284 ) {
1285 mesh_boolean_operation(result, A, B, "A-B", verbose);
1286 }
1287
1295 Mesh& M, index_t max_iter = 3, bool verbose=false
1296 );
1297
1308 Mesh& M, index_t f1, index_t f2
1309 );
1310
1311 /**************************************************************************/
1312}
1313
1314#endif
#define geo_assert_not_reached
Sets a non reachable point in the program.
Definition assert.h:177
#define geo_assert(x)
Verifies that a condition is met.
Definition assert.h:149
#define geo_debug_assert(x)
Verifies that a condition is met.
Definition assert.h:196
Generic mechanism for attributes.
Manages an attribute attached to a set of object.
Detects and retriangulates a set of coplanar facets for MeshSurfaceIntersection.
index_t adjacent_facet(index_t c) const
Gets the facet that a corner is adjacent to.
Definition mesh.h:951
void set_adjacent_facet(index_t c, index_t f)
Sets the facet that a corner is adjacent to.
Definition mesh.h:1014
index_t vertex(index_t f, index_t lv) const
Gets a vertex by facet and local vertex index.
Definition mesh.h:1143
index_range corners(index_t f) const
Gets the corners of a facet.
Definition mesh.h:1538
Meshes a single triangle with the constraints that come from the intersections with the other triangl...
AttributesManager & attributes() const
Gets the attributes manager.
Definition mesh.h:109
index_t nb() const
Gets the number of (sub-)elements.
Definition mesh.h:98
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.
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.
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_fine_verbose(bool x)
Display detailed information while computing the intersection. Default is unset.
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...
Represents a mesh.
Definition mesh.h:3068
Wraps an integer to be used with the range-based for construct.
Definition range.h:66
Comparator class for vec3Hg \detail Used to create maps indexed by vec3Hg or SOS symbolic perturbatio...
Definition vechg.h:293
3d vector with homogeneous coordinates
Definition vechg.h:185
Vector with aligned memory allocation.
Definition memory.h:660
index_t size() const
Gets the number of elements.
Definition memory.h:706
Exact predicates and constructs.
Common include file, providing basic definitions. Should be included before anything else by all head...
The class that represents a mesh.
Functions to load and save meshes.
SOSMode
Mode for symbolic perturbations.
Definition predicates.h:67
std::atomic_flag spinlock
A lightweight synchronization structure.
Global Vorpaline namespace.
Definition basic.h:55
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_intersection(Mesh &result, const Mesh &A, const Mesh &B, MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT)
Computes the intersection of two surface meshes.
void mesh_difference(Mesh &result, const Mesh &A, const Mesh &B, MeshBooleanOperationFlags flags=MESH_BOOL_OPS_DEFAULT)
Computes the difference of two surface meshes.
void mesh_boolean_operation(Mesh &result, const Mesh &A, const 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, const Mesh &A, const 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.
Definition numeric.h:329
Sign
Integer constants that represent the sign of a value.
Definition numeric.h:68
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...
Function and classes for process manipulation.