Iterative epsilon in faceset helper

src/ifcgeom/IfcGeomFunctions.cpp

@@ -4137,6 +4137,14 @@ IfcGeom::Kernel::faceset_helper::faceset_helper(Kernel* kernel, const IfcSchema:
 	// Use the bbox diagonal to influence local epsilon
 	// double bdiff = std::sqrt(box.SquareExtent());
 
+	// @todo the bounding box diagonal is not used (see above)
+	// because we're explicitly interested in the miminal
+	// dimension of the element to limit the tolerance (for sheet-
+	// like elements for example). But the way below is very
+	// dependent on orientation due to the usage of the
+	// axis-aligned bounding box. Use PCA to find three non-aligned
+	// set of dimensions and use the one with the smallest eigenvalue.
+
 	// Find the minimal bounding box edge
 	double bmin[3], bmax[3];
 	box.Get(bmin[0], bmin[1], bmin[2], bmax[0], bmax[1], bmax[2]);
@@ -4150,15 +4158,37 @@ IfcGeom::Kernel::faceset_helper::faceset_helper(Kernel* kernel, const IfcSchema:
 
 	eps_ = kernel->getValue(GV_PRECISION) * 10. * (std::min)(1.0, bdiff);
 
+	// @todo, there a tiny possibility that the duplicate faces are triggered
+	// for an internal boundary, that is also present as an external boundary.
+	// This will result in non-manifold configuration then, but this is deemed
+	// such as corner-case that it is not considered.
+	IfcSchema::IfcPolyLoop::list::ptr loops = IfcParse::traverse((IfcUtil::IfcBaseClass*)l)->as<IfcSchema::IfcPolyLoop>();
+
+	size_t loops_removed, non_manifold, duplicate_faces;
+
+	std::map<std::pair<int, int>, int> edge_use;
+
+	for (int i = 0; i < 3; ++i) {
+		// Some times files, have large tolerance values specified collapsing too many vertices.
+		// This case we detect below and re-run the loop with smaller epsilon. Normally
+		// the body of this loop would only be executed once.
+
+		loops_removed = 0;
+		non_manifold = 0;
+		duplicate_faces = 0;
+
+		vertex_mapping_.clear();
+		duplicates_.clear();
+
+		edge_use.clear();
+
 		if (eps_ < Precision::Confusion()) {
 			// occt uses some hard coded precision values, don't go smaller than that.
 			// @todo, can be reset though with BRepLib::Precision(double)
 			eps_ = Precision::Confusion();
 		}		
 
-	std::map<std::pair<int, int>, int> edge_use;
+		for (int i = 0; i < (int)pnts.size(); ++i) {
-
-	for (int i = 0; i < (int) pnts.size(); ++i) {
 			if (pnts[i]) {
 				std::set<int> vs;
 				find_neighbours(tree, pnts, vs, i, eps_);
@@ -4171,14 +4201,6 @@ IfcGeom::Kernel::faceset_helper::faceset_helper(Kernel* kernel, const IfcSchema:
 			}
 		}		
 
-	// @todo, there a tiny possibility that the duplicate faces are triggered
-	// for an internal boundary, that is also present as an external boundary.
-	// This will result in non-manifold configuration then, but this is deemed
-	// such as corner-case that it is not considered.
-	IfcSchema::IfcPolyLoop::list::ptr loops = IfcParse::traverse((IfcUtil::IfcBaseClass*)l)->as<IfcSchema::IfcPolyLoop>();
-
-	size_t loops_removed = 0, non_manifold = 0, duplicate_faces = 0;
-
 		typedef std::array<int, 2> edge_t;
 		typedef std::set<edge_t> edge_set_t;
 		std::set<edge_set_t> edge_sets;
@@ -4190,7 +4212,7 @@ IfcGeom::Kernel::faceset_helper::faceset_helper(Kernel* kernel, const IfcSchema:
 			edge_set_t segment_set;
 
 			loop_(ps, [&segments, &segment_set](int C, int D, bool) {
-			segment_set.insert({{ C, D }});
+				segment_set.insert({ { C, D } });
 				segments.push_back({ C, D });
 			});
 
@@ -4210,6 +4232,13 @@ IfcGeom::Kernel::faceset_helper::faceset_helper(Kernel* kernel, const IfcSchema:
 			}
 		}
 
+		if (edge_use.size() != 0) {
+			break;
+		} else {
+			eps_ /= 10.;
+		}
+	}
+
 	for (auto& p : edge_use) {
 		int a, b;
 		std::tie(a, b) = p.first;