Missing files

src/ifcgeom_schema_agnostic/IteratorCache.h

@@ -0,0 +1,31 @@
+/********************************************************************************
+ *                                                                              *
+ * This file is part of IfcOpenShell.                                           *
+ *                                                                              *
+ * IfcOpenShell is free software: you can redistribute it and/or modify         *
+ * it under the terms of the Lesser GNU General Public License as published by  *
+ * the Free Software Foundation, either version 3.0 of the License, or          *
+ * (at your option) any later version.                                          *
+ *                                                                              *
+ * IfcOpenShell is distributed in the hope that it will be useful,              *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of               *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the                 *
+ * Lesser GNU General Public License for more details.                          *
+ *                                                                              *
+ * You should have received a copy of the Lesser GNU General Public License     *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>.         *
+ *                                                                              *
+ ********************************************************************************/
+
+#ifndef ITERATOR_CACHE_H
+#define ITERATOR_CACHE_H
+
+namespace IfcGeom {
+
+	class IteratorCache {
+
+	};
+
+}
+
+#endif

src/ifcgeom_schema_agnostic/boolean_utils.cpp

@@ -0,0 +1,777 @@
+#include "boolean_utils.h"
+
+#include "../ifcgeom/IfcGeomTree.h"
+
+#include <BRepBuilderAPI_Copy.hxx>
+#include <TopExp_Explorer.hxx>
+#include <GProp_GProps.hxx>
+#include <BRepGProp.hxx>
+#include <TopExp.hxx>
+#include <TopoDS.hxx>
+#include <Bnd_Box.hxx>
+#include <Extrema_ExtPC.hxx>
+#include <Geom_Plane.hxx>
+#include <Geom_BSplineCurve.hxx>
+#include <ShapeUpgrade_UnifySameDomain.hxx>
+#include <GeomAPI_ExtremaCurveCurve.hxx>
+#include <ShapeAnalysis_Surface.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+
+#include <vector>
+
+void IfcGeom::util::copy_operand(const TopTools_ListOfShape & l, TopTools_ListOfShape & r) {
+#if OCC_VERSION_HEX < 0x70000
+	TopTools_ListIteratorOfListOfShape it(l);
+	for (; it.More(); it.Next()) {
+		r.Append(BRepBuilderAPI_Copy(it.Value()));
+	}
+#else
+	// On OCCT 7.0 and higher BRepAlgoAPI_BuilderAlgo::SetNonDestructive(true) is
+	// called. Not entirely sure on the behaviour before 7.0, so overcautiously
+	// create copies.
+	r.Assign(l);
+#endif
+}
+
+TopoDS_Shape IfcGeom::util::copy_operand(const TopoDS_Shape & s) {
+#if OCC_VERSION_HEX < 0x70000
+	return BRepBuilderAPI_Copy(s);
+#else
+	return s;
+#endif
+}
+
+double IfcGeom::util::min_edge_length(const TopoDS_Shape & a) {
+	double min_edge_len = std::numeric_limits<double>::infinity();
+	TopExp_Explorer exp(a, TopAbs_EDGE);
+	for (; exp.More(); exp.Next()) {
+		GProp_GProps prop;
+		BRepGProp::LinearProperties(exp.Current(), prop);
+		double l = prop.Mass();
+		if (l < min_edge_len) {
+			min_edge_len = l;
+		}
+	}
+	return min_edge_len;
+}
+
+double IfcGeom::util::min_vertex_edge_distance(const TopoDS_Shape & a, double min_search, double max_search) {
+	double M = std::numeric_limits<double>::infinity();
+
+	TopTools_IndexedMapOfShape vertices, edges;
+
+	TopExp::MapShapes(a, TopAbs_VERTEX, vertices);
+	TopExp::MapShapes(a, TopAbs_EDGE, edges);
+
+	IfcGeom::impl::tree<int> tree;
+
+	// Add edges to tree
+	for (int i = 1; i <= edges.Extent(); ++i) {
+		tree.add(i, edges(i));
+	}
+
+	for (int j = 1; j <= vertices.Extent(); ++j) {
+		const TopoDS_Vertex& v = TopoDS::Vertex(vertices(j));
+		gp_Pnt p = BRep_Tool::Pnt(v);
+
+		Bnd_Box b;
+		b.Add(p);
+		b.Enlarge(max_search);
+
+		std::vector<int> edge_idxs = tree.select_box(b, false);
+		std::vector<int>::const_iterator it = edge_idxs.begin();
+		for (; it != edge_idxs.end(); ++it) {
+			const TopoDS_Edge& e = TopoDS::Edge(edges(*it));
+			TopoDS_Vertex v1, v2;
+			TopExp::Vertices(e, v1, v2);
+
+			if (v.IsSame(v1) || v.IsSame(v2)) {
+				continue;
+			}
+
+			BRepAdaptor_Curve crv(e);
+			Extrema_ExtPC ext(p, crv);
+			if (!ext.IsDone()) {
+				continue;
+			}
+
+			for (int i = 1; i <= ext.NbExt(); ++i) {
+				const double m = sqrt(ext.SquareDistance(i));
+				if (m < M && m > min_search) {
+					M = m;
+				}
+			}
+		}
+	}
+
+	return M;
+}
+
+bool IfcGeom::util::faces_overlap(const TopoDS_Face & f, const TopoDS_Face & g) {
+	points_on_planar_face_generator pgen(f);
+
+	BRep_Builder B;
+	gp_Pnt test;
+	double eps = BRep_Tool::Tolerance(f) + BRep_Tool::Tolerance(g);
+
+	BRepExtrema_DistShapeShape x;
+	x.LoadS1(g);
+
+	while (pgen(test)) {
+		TopoDS_Vertex V;
+		B.MakeVertex(V, test, Precision::Confusion());
+		x.LoadS2(V);
+		x.Perform();
+		if (x.IsDone() && x.NbSolution() == 1) {
+			if (x.Value() > eps) {
+				return false;
+			}
+		}
+	}
+
+	return true;
+}
+
+double IfcGeom::util::min_face_face_distance(const TopoDS_Shape & a, double max_search) {
+	/*
+	NB: This is currently only implemented for planar surfaces.
+	*/
+	double M = std::numeric_limits<double>::infinity();
+
+	TopTools_IndexedMapOfShape faces;
+
+	TopExp::MapShapes(a, TopAbs_FACE, faces);
+
+	IfcGeom::impl::tree<int> tree;
+
+	// Add faces to tree
+	for (int i = 1; i <= faces.Extent(); ++i) {
+		if (BRep_Tool::Surface(TopoDS::Face(faces(i)))->DynamicType() == STANDARD_TYPE(Geom_Plane)) {
+			tree.add(i, faces(i));
+		}
+	}
+
+	for (int j = 1; j <= faces.Extent(); ++j) {
+		const TopoDS_Face& f = TopoDS::Face(faces(j));
+		const Handle(Geom_Surface)& fs = BRep_Tool::Surface(f);
+
+		if (fs->DynamicType() != STANDARD_TYPE(Geom_Plane)) {
+			continue;
+		}
+
+		points_on_planar_face_generator pgen(f);
+
+		Bnd_Box b;
+		BRepBndLib::AddClose(f, b);
+		b.Enlarge(max_search);
+
+		std::vector<int> face_idxs = tree.select_box(b, false);
+		std::vector<int>::const_iterator it = face_idxs.begin();
+		for (; it != face_idxs.end(); ++it) {
+			if (*it == j) {
+				continue;
+			}
+
+			const TopoDS_Face& g = TopoDS::Face(faces(*it));
+			const Handle(Geom_Surface)& gs = BRep_Tool::Surface(g);
+
+			auto p0 = Handle(Geom_Plane)::DownCast(fs);
+			auto p1 = Handle(Geom_Plane)::DownCast(gs);
+
+			if (p0->Position().IsCoplanar(p1->Position(), max_search, asin(max_search))) {
+				pgen.reset();
+
+				BRepTopAdaptor_FClass2d cls(g, BRep_Tool::Tolerance(g));
+
+				gp_Pnt test;
+				while (pgen(test)) {
+					gp_Vec d = test.XYZ() - p1->Position().Location().XYZ();
+					double u = d.Dot(p1->Position().XDirection());
+					double v = d.Dot(p1->Position().YDirection());
+
+					// nb: TopAbs_ON is explicitly not considered to prevent matching adjacent faces
+					// with similar orientations.
+					if (cls.Perform(gp_Pnt2d(u, v)) == TopAbs_IN) {
+						gp_Pnt test2;
+						p1->D0(u, v, test2);
+						double w = std::abs(gp_Vec(p1->Position().Direction().XYZ()).Dot(test2.XYZ() - test.XYZ()));
+						if (w < M) {
+							M = w;
+						}
+					}
+				}
+			}
+		}
+	}
+
+	return M;
+}
+
+int IfcGeom::util::bounding_box_overlap(double p, const TopoDS_Shape & a, const TopTools_ListOfShape & b, TopTools_ListOfShape & c) {
+	int N = 0;
+
+	Bnd_Box A;
+	BRepBndLib::Add(a, A);
+
+	if (A.IsVoid()) {
+		return 0;
+	}
+
+	TopTools_ListIteratorOfListOfShape it(b);
+	for (; it.More(); it.Next()) {
+		Bnd_Box B;
+		BRepBndLib::Add(it.Value(), B);
+
+		if (B.IsVoid()) {
+			continue;
+		}
+
+		if (A.Distance(B) < p) {
+			c.Append(it.Value());
+		} else {
+			++N;
+		}
+	}
+
+	return N;
+}
+
+bool IfcGeom::util::get_edge_axis(const TopoDS_Edge & e, gp_Ax1 & ax) {
+	double _, __;
+
+	auto crv = BRep_Tool::Curve(e, _, __);
+	auto line = Handle_Geom_Line::DownCast(crv);
+	auto bsple = Handle_Geom_BSplineCurve::DownCast(crv);
+
+	if (line) {
+		ax = line->Position();
+		return true;
+	} else if (bsple) {
+		if (bsple->NbPoles() == 2 && bsple->Degree() == 1) {
+			gp_Dir V(bsple->Poles().Last().XYZ() - bsple->Poles().First().XYZ());
+			ax = gp_Ax1(bsple->Poles().First(), V);
+			return true;
+		}
+	}
+
+	return false;
+}
+
+bool IfcGeom::util::is_subset(const TopTools_IndexedMapOfShape & lhs, const TopTools_IndexedMapOfShape & rhs) {
+	if (rhs.Extent() < lhs.Extent()) {
+		return false;
+	}
+	for (int i = 1; i < lhs.Extent(); ++i) {
+		auto& s = lhs.FindKey(i);
+		if (!rhs.Contains(s)) {
+			return false;
+		}
+	}
+	return true;
+}
+
+bool IfcGeom::util::is_extrusion(const gp_Vec & v, const TopoDS_Shape & s, TopoDS_Face & base, std::pair<double, double>& interval) {
+	// This assumes UnifySameDomain has been processed on s, so that
+	// the extrusion top and bottom are a single face.
+
+	TopTools_IndexedDataMapOfShapeListOfShape mapping;
+	TopExp::MapShapesAndAncestors(s, TopAbs_EDGE, TopAbs_FACE, mapping);
+	TopExp::MapShapesAndAncestors(s, TopAbs_VERTEX, TopAbs_FACE, mapping);
+
+	TopTools_ListOfShape parallel;
+	TopTools_IndexedMapOfShape curved_orthogonal;
+	gp_Ax1 ax;
+	gp_Ax1 V(gp::Origin(), v);
+
+	// Segment edges in parallel to extrusion direction, and orthogonal or curved,
+	// where the latter two categories have to make the edges part of the base or
+	// top face. When neither of these categories the shape is not a extrusion
+	// or the extrusion direction is not orthogonal to its basis.
+	for (int i = 1; i < mapping.Extent(); ++i) {
+		auto& s = mapping.FindKey(i);
+		if (s.ShapeType() != TopAbs_EDGE) {
+			continue;
+		}
+
+		// @todo use a linear tolernace and the face extrimities, see #2218
+		const TopoDS_Edge& e = TopoDS::Edge(s);
+		if (!get_edge_axis(e, ax)) {
+			// curved
+			curved_orthogonal.Add(e);
+		} else if (ax.IsParallel(V, 1.e-7)) {
+			parallel.Append(e);
+		} else if (ax.IsNormal(V, 1.e-7)) {
+			// ortho
+			curved_orthogonal.Add(e);
+		} else {
+			return false;
+		}
+	}
+
+	// Select the two faces for which their edges are subsets
+	// of the ortho/curved edges
+	TopTools_IndexedMapOfShape ortho_faces;
+	for (TopExp_Explorer exp(s, TopAbs_FACE); exp.More(); exp.Next()) {
+		TopTools_IndexedMapOfShape face_edges;
+		TopExp::MapShapes(exp.Current(), TopAbs_EDGE, face_edges);
+		if (is_subset(face_edges, curved_orthogonal)) {
+			ortho_faces.Add(exp.Current());
+		}
+	}
+
+	// There should be a basis and top face
+	if (ortho_faces.Extent() != 2) {
+		return false;
+	}
+
+	// For the parallel edges assert that its two vertices are part
+	// of both the basis and the top face.
+	for (TopTools_ListIteratorOfListOfShape it(parallel);
+		it.More(); it.Next()) {
+		TopoDS_Vertex v01[2];
+		TopExp::Vertices(TopoDS::Edge(it.Value()), v01[0], v01[1]);
+
+		TopTools_IndexedMapOfShape v_ortho_faces;
+		int nb_ortho_faces[2] = { 0,0 };
+
+		for (int i = 0; i < 2; ++i) {
+			auto& faces = mapping.FindFromKey(v01[i]);
+
+			for (TopTools_ListIteratorOfListOfShape jt(faces);
+				jt.More(); jt.Next()) {
+				if (ortho_faces.Contains(jt.Value())) {
+					nb_ortho_faces[i] ++;
+					v_ortho_faces.Add(jt.Value());
+				}
+			}
+		}
+
+		bool sets_equal = v_ortho_faces.Size() == ortho_faces.Size() && is_subset(v_ortho_faces, ortho_faces);
+		if (!sets_equal) {
+			return false;
+		}
+	}
+
+	// Assert the base/top faces are planar and get the interval
+	// (dot products along axis) for which the extrusion is defined
+	// If necessary swap the two faces so that the basis face has
+	// the smallest dot product along the axis.
+	auto f0 = TopoDS::Face(ortho_faces.FindKey(1));
+	auto f1 = TopoDS::Face(ortho_faces.FindKey(2));
+
+	const Handle(Geom_Surface)& f0_s = BRep_Tool::Surface(f0);
+	const Handle(Geom_Surface)& f1_s = BRep_Tool::Surface(f1);
+
+	auto p0 = Handle(Geom_Plane)::DownCast(f0_s);
+	auto p1 = Handle(Geom_Plane)::DownCast(f1_s);
+
+	if (p0.IsNull() || p1.IsNull()) {
+		return false;
+	}
+
+	auto dot0 = p0->Location().XYZ().Dot(v.XYZ());
+	auto dot1 = p1->Location().XYZ().Dot(v.XYZ());
+
+	if (dot0 > dot1) {
+		std::swap(dot0, dot1);
+		std::swap(f0, f1);
+	}
+
+	base = f0;
+	interval = { dot0, dot1 };
+
+	return true;
+}
+
+int IfcGeom::util::eliminate_touching_operands(double prec, const TopoDS_Shape & a, const TopTools_ListOfShape & bs, TopTools_ListOfShape & c) {
+	TopTools_IndexedMapOfShape a_faces;
+	TopExp::MapShapes(a, TopAbs_FACE, a_faces);
+
+	// Check if any of the faces in a are non-planar, which is
+	// not supported by this quick check.
+	for (int i = 1; i <= a_faces.Extent(); ++i) {
+		auto surf = BRep_Tool::Surface(TopoDS::Face(a_faces(i)));
+		if (surf->DynamicType() != STANDARD_TYPE(Geom_Plane)) {
+			return 0;
+		}
+	}
+
+	TopTools_IndexedMapOfShape a_vertices;
+	TopExp::MapShapes(a, TopAbs_VERTEX, a_vertices);
+
+	IfcGeom::impl::tree<int> tree;
+
+	// Add faces to tree
+	for (int i = 1; i <= a_faces.Extent(); ++i) {
+		tree.add(i, a_faces(i));
+	}
+
+	int N = 0;
+
+	TopTools_ListIteratorOfListOfShape it(bs);
+	for (; it.More(); it.Next()) {
+		bool is_touching = false;
+
+		auto& b = it.Value();
+
+		TopTools_IndexedMapOfShape b_faces;
+		TopExp::MapShapes(b, TopAbs_FACE, b_faces);
+
+		// Check if any of the faces in b are non-planar, which is
+		// not supported by this quick check.
+		for (int i = 1; i <= b_faces.Extent(); ++i) {
+			auto surf = BRep_Tool::Surface(TopoDS::Face(b_faces(i)));
+			if (surf->DynamicType() != STANDARD_TYPE(Geom_Plane)) {
+				continue;
+			}
+		}
+
+		TopTools_IndexedMapOfShape b_vertices;
+		TopExp::MapShapes(b, TopAbs_VERTEX, b_vertices);
+
+		for (int k = 1; k <= b_faces.Extent(); ++k) {
+			const TopoDS_Face& f_b = TopoDS::Face(b_faces(k));
+			Bnd_Box B;
+			BRepBndLib::Add(f_b, B);
+
+			// Query tree using b_face bounding box
+			for (auto& i : tree.select_box(B, false)) {
+				const TopoDS_Face& f_a = TopoDS::Face(a_faces(i));
+
+				TopTools_IndexedMapOfShape f_a_vertices;
+				TopExp::MapShapes(f_a, TopAbs_VERTEX, f_a_vertices);
+
+				BRepGProp_Face prop_a(f_a);
+				BRepGProp_Face prop_b(f_b);
+
+				gp_Pnt p_a, p_b;
+				gp_Vec v_a, v_b;
+
+				double u0, u1, v0, v1;
+				prop_a.Bounds(u0, u1, v0, v1);
+				prop_a.Normal((u0 + u1) / 2., (u0 + u1) / 2., p_a, v_a);
+
+				prop_b.Bounds(u0, u1, v0, v1);
+				prop_b.Normal((u0 + u1) / 2., (u0 + u1) / 2., p_b, v_b);
+
+				bool all_vertices_behind_f_a = true;
+
+				// Check if all 'other' vertices in a are pointing
+				// away from the face in a, so that there is no geometry
+				// from a in front of the face that could participate
+				// in the boolean subtraction.
+				for (int j = 1; j <= a_vertices.Extent(); ++j) {
+					if (!f_a_vertices.Contains(a_vertices(j))) {
+						auto p = BRep_Tool::Pnt(TopoDS::Vertex(a_vertices(j)));
+						if ((p.XYZ() - p_a.XYZ()).Dot(v_a.XYZ()) > prec) {
+							all_vertices_behind_f_a = false;
+							break;
+						}
+					}
+				}
+
+				if (!all_vertices_behind_f_a) {
+					continue;
+				}
+
+				// Check if surface normals are opposite
+				if (v_a.IsOpposite(v_b, 1.e-5)) {
+					// Check if faces are co-planar
+					if ((p_b.XYZ() - p_a.XYZ()).Dot(v_a.XYZ()) <= prec) {
+
+						TopTools_IndexedMapOfShape f_b_vertices;
+						TopExp::MapShapes(f_b, TopAbs_VERTEX, f_b_vertices);
+
+						bool all_vertices_behind_f_b = true;
+
+						// Check if all 'other' vertices in b are pointing
+						// away from the face in a. So that a boolean subtraction
+						// would not alter a.
+						for (int j = 1; j <= b_vertices.Extent(); ++j) {
+							if (!f_b_vertices.Contains(b_vertices(j))) {
+								auto p = BRep_Tool::Pnt(TopoDS::Vertex(b_vertices(j)));
+								if ((p.XYZ() - p_a.XYZ()).Dot(v_a.XYZ()) < prec * 10.) {
+									all_vertices_behind_f_b = false;
+									break;
+								}
+							}
+						}
+
+						if (all_vertices_behind_f_b) {
+							is_touching = true;
+							break;
+						}
+
+					}
+				}
+			}
+
+			if (is_touching) {
+				break;
+			}
+		}
+
+		if (!is_touching) {
+			c.Append(it.Value());
+		} else {
+			++N;
+		}
+	}
+
+	return N;
+}
+
+TopoDS_Shape IfcGeom::util::unify(const TopoDS_Shape & s, double tolerance) {
+	tolerance = (std::min)(min_edge_length(s) / 2., tolerance);
+	ShapeUpgrade_UnifySameDomain usd(s);
+#if OCC_VERSION_HEX >= 0x70200
+	usd.SetSafeInputMode(true);
+#endif
+#if OCC_VERSION_HEX >= 0x70100
+	usd.SetLinearTolerance(tolerance);
+	usd.SetAngularTolerance(1.e-3);
+#endif
+	usd.Build();
+	return usd.Shape();
+}
+
+bool IfcGeom::util::boolean_subtraction_2d_using_builder(const TopoDS_Shape & a_input, const TopTools_ListOfShape & b_input, TopoDS_Shape & result, double eps) {
+	IfcGeom::impl::tree<int> edge_tree;
+
+	TopTools_ListOfShape ab_input = b_input;
+	ab_input.Prepend(a_input);
+
+	TopTools_ListIteratorOfListOfShape it(ab_input);
+	int shape_index = 0;
+	int edge_index = 0;
+	std::map<int, int> edge_index_to_shape_index;
+
+	std::vector<TopoDS_Shape> shapes;
+	std::vector<std::pair<size_t, TopoDS_Edge>> edges;
+	// First is the outer wire
+	std::vector<TopoDS_Wire> wires;
+
+	for (; it.More(); it.Next(), ++shape_index) {
+		if (it.Value().ShapeType() != TopAbs_FACE) {
+			return false;
+		}
+
+		const TopoDS_Face& f = TopoDS::Face(it.Value());
+		TopoDS_Wire outer_wire;
+
+		if (shape_index == 0) {
+			outer_wire = BRepTools::OuterWire(f);
+			wires.push_back(outer_wire);
+		}
+
+		size_t num_wires = 0;
+		TopoDS_Iterator it2(it.Value());
+		for (; it2.More(); it2.Next()) {
+			++num_wires;
+
+			if (outer_wire.IsNull() || !it2.Value().IsSame(outer_wire)) {
+				wires.push_back(TopoDS::Wire(it2.Value()));
+
+				if (shape_index == 0 && num_wires > 0) {
+					// An inner wire on the first operand face: reverse, because
+					// MakeFace expects inner boundaries to be added as bounded
+					// areas.
+					wires.back().Reverse();
+				}
+			}
+		}
+
+		if (num_wires > 1 && shape_index != 0) {
+			// The first operand can have inner wires, but the others
+			// can't because a inner wire would result in an additional
+			// outer wire for the result.
+			return false;
+		}
+
+		shapes.push_back(it.Value());
+		TopExp_Explorer exp(it.Value(), TopAbs_EDGE);
+		for (; exp.More(); exp.Next(), ++edge_index) {
+			edge_tree.add(edge_index, exp.Current());
+			edge_index_to_shape_index[edge_index] = shape_index;
+			edges.push_back({ shape_index, TopoDS::Edge(exp.Current()) });
+		}
+	}
+
+	shape_index = 0;
+	edge_index = 0;
+
+	it.Initialize(ab_input);
+	for (; it.More(); it.Next(), ++shape_index) {
+		TopExp_Explorer exp(it.Value(), TopAbs_EDGE);
+		for (; exp.More(); exp.Next(), ++edge_index) {
+			Bnd_Box b;
+			BRepBndLib::Add(exp.Current(), b);
+			b.Enlarge(eps);
+
+			for (auto& i : edge_tree.select_box(b)) {
+				if (i == edge_index) {
+					// Skip self-selection
+					continue;
+				}
+
+				if (edges[i].first == shape_index) {
+					// Skip edges of the same operand
+					continue;
+				}
+
+				const TopoDS_Edge& e0 = TopoDS::Edge(exp.Current());
+				const TopoDS_Edge& e1 = edges[i].second;
+
+				double u11, u12, u21, u22, U1, U2;
+
+				GeomAPI_ExtremaCurveCurve ecc(
+					BRep_Tool::Curve(e0, u11, u12),
+					BRep_Tool::Curve(e1, u21, u22)
+				);
+
+				// @todo: extend this to work in case of multiple extrema and curved segments.
+				const bool unbounded_intersects = (!ecc.Extrema().IsParallel() && ecc.NbExtrema() == 1 && ecc.Distance(1) < eps);
+				if (unbounded_intersects) {
+					ecc.Parameters(1, U1, U2);
+
+					if (u11 > u12) {
+						std::swap(u11, u12);
+					}
+					if (u21 > u22) {
+						std::swap(u21, u22);
+					}
+
+					/// @todo: tfk: probably need different thresholds on non-linear curves
+					u11 -= eps;
+					u12 += eps;
+					u21 -= eps;
+					u22 += eps;
+
+					if (u11 < U1 && U1 < u12 && u21 < U2 && U2 < u22) {
+						// Edge curves belonging to different operands intersect, don't process
+						// using builder.
+						Logger::Notice("Intersecting boundaries");
+						return false;
+					}
+				}
+			}
+		}
+	}
+
+	// Only inner wires are considered that are directly contained in the outer wire
+	// Redundant subtractions are eliminated.
+
+	std::vector<bool> redundant(wires.size(), false);
+
+	std::vector<TopoDS_Face> wire_faces;
+	wire_faces.reserve(wires.size());
+
+	std::vector<BRepTopAdaptor_FClass2d> wire_clss;
+	wire_clss.reserve(wires.size());
+
+	std::vector<std::unique_ptr<ShapeAnalysis_Surface>> sass;
+	sass.reserve(wires.size());
+
+	for (auto& w : wires) {
+		wire_faces.push_back(BRepBuilderAPI_MakeFace(w).Face());
+		wire_clss.emplace_back(wire_faces.back(), eps);
+		sass.push_back(std::make_unique<ShapeAnalysis_Surface>(BRep_Tool::Surface(wire_faces.back())));
+	}
+
+	// First check for containment in outer wire
+	for (auto it = ++wires.begin(); it != wires.end(); ++it) {
+		// Considering a single vertex is sufficient because we have already
+		// guaranteed that the edges of different operands do not cross.
+		TopoDS_Iterator it_ed(*it);
+		auto& ed = it_ed.Value();
+
+		TopoDS_Iterator it_v(ed);
+		auto& v = TopoDS::Vertex(it_v.Value());
+
+		auto pnt = BRep_Tool::Pnt(v);
+		auto p2d = sass[0]->ValueOfUV(pnt, eps);
+		if (wire_clss[0].Perform(p2d) != TopAbs_IN) {
+			// A wire is not contained in the outer wire, it's a subtraction without
+			// any effect and marked as redundant. Feeding it to the builder algo
+			// will likely cause problems.
+			redundant[std::distance(wires.begin(), it)] = true;
+			Logger::Notice("Subtraction operand outside of outer bound");
+		}
+	}
+
+	// Now build a tree to find inner wires contained in other inner wires
+	// NB first wire is *not* in this tree
+	IfcGeom::impl::tree<int> wire_tree;
+	for (size_t wire_index = 1; wire_index < wires.size(); ++wire_index) {
+		wire_tree.add(wire_index, wires[wire_index]);
+	}
+
+	for (size_t wire_index = 1; wire_index < wires.size(); ++wire_index) {
+		Bnd_Box b;
+		BRepBndLib::Add(wires[wire_index], b);
+		b.Enlarge(eps);
+
+		// We're only selecting operands completely within b because we
+		// have already guaranteed they do not intersect. So they are
+		// either fully in or out. Selecting with complete_within=true
+		// will filter out some unnecessary cases. It also means we need
+		// that due this assymetry we need to process all pairs of wire
+		// indices and not just the pairs where the first element is less
+		// than the second element.
+		for (auto& other_index : wire_tree.select_box(b, true)) {
+			// other_index is fully contained in wire_index
+			if (wire_index == other_index) {
+				continue;
+			}
+
+			TopoDS_Iterator it_ed(wires[other_index]);
+			auto& ed = it_ed.Value();
+
+			TopoDS_Iterator it_v(ed);
+			auto& v = TopoDS::Vertex(it_v.Value());
+
+			auto pnt = BRep_Tool::Pnt(v);
+			auto p2d = sass[wire_index]->ValueOfUV(pnt, eps);
+			if (wire_clss[wire_index].Perform(p2d) == TopAbs_IN) {
+				// A wire is contained within another operand
+				redundant[other_index] = true;
+				Logger::Notice("Subtraction operand contained in other");
+			}
+		}
+	}
+
+	BRepBuilderAPI_MakeFace mf(wire_faces[0]);
+	for (size_t wire_index = 1; wire_index < wires.size(); ++wire_index) {
+		if (!redundant[wire_index]) {
+			mf.Add(TopoDS::Wire(wires[wire_index].Reversed()));
+		}
+	}
+	result = mf.Face();
+
+	return true;
+}
+
+void IfcGeom::util::points_on_planar_face_generator::reset() {
+	i = j = (int)inset_;
+}
+
+bool IfcGeom::util::points_on_planar_face_generator::operator()(gp_Pnt& p) {
+	while (j < N) {
+		double u = u0 + (u1 - u0) * i / N;
+		double v = v0 + (v1 - v0) * j / N;
+
+		i++;
+		if (i == N) {
+			i = 0;
+			j++;
+		}
+
+		// Specifically does not consider ON
+		if (cls_.Perform(gp_Pnt2d(u, v)) == TopAbs_IN) {
+			plane_->D0(u, v, p);
+			return true;
+		}
+	}
+
+	return false;
+}

src/ifcgeom_schema_agnostic/boolean_utils.h

@@ -0,0 +1,90 @@
+/********************************************************************************
+ *                                                                              *
+ * This file is part of IfcOpenShell.                                           *
+ *                                                                              *
+ * IfcOpenShell is free software: you can redistribute it and/or modify         *
+ * it under the terms of the Lesser GNU General Public License as published by  *
+ * the Free Software Foundation, either version 3.0 of the License, or          *
+ * (at your option) any later version.                                          *
+ *                                                                              *
+ * IfcOpenShell is distributed in the hope that it will be useful,              *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of               *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the                 *
+ * Lesser GNU General Public License for more details.                          *
+ *                                                                              *
+ * You should have received a copy of the Lesser GNU General Public License     *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>.         *
+ *                                                                              *
+ ********************************************************************************/
+
+#ifndef BOOLEAN_UTILS_H
+#define BOOLEAN_UTILS_H
+
+#include <TopoDS_Shape.hxx>
+#include <TopTools_ListOfShape.hxx>
+#include <Geom_Surface.hxx>
+#include <TopoDS_Face.hxx>
+#include <BRepTopAdaptor_FClass2d.hxx>
+#include <BRep_Tool.hxx>
+#include <BRepTools.hxx>
+#include <TopTools_IndexedMapOfShape.hxx>
+
+namespace IfcGeom {
+	namespace util {
+
+		void copy_operand(const TopTools_ListOfShape& l, TopTools_ListOfShape& r);
+
+		TopoDS_Shape copy_operand(const TopoDS_Shape& s);
+
+		double min_edge_length(const TopoDS_Shape& a);
+
+		double min_vertex_edge_distance(const TopoDS_Shape& a, double min_search, double max_search);
+
+		class points_on_planar_face_generator {
+		private:
+			const TopoDS_Face& f_;
+			Handle(Geom_Surface) plane_;
+			BRepTopAdaptor_FClass2d cls_;
+			double u0, u1, v0, v1;
+			int i, j;
+			bool inset_;
+			static const int N = 10;
+
+		public:
+			points_on_planar_face_generator(const TopoDS_Face& f, bool inset = false)
+				: f_(f)
+				, plane_(BRep_Tool::Surface(f_))
+				, cls_(f_, BRep_Tool::Tolerance(f_))
+				, i((int)inset), j((int)inset)
+				, inset_(inset)
+			{
+				BRepTools::UVBounds(f_, u0, u1, v0, v1);
+			}
+
+			void reset();
+
+			bool operator()(gp_Pnt& p);
+		};
+
+		bool faces_overlap(const TopoDS_Face& f, const TopoDS_Face& g);
+
+		double min_face_face_distance(const TopoDS_Shape& a, double max_search);
+
+		int bounding_box_overlap(double p, const TopoDS_Shape& a, const TopTools_ListOfShape& b, TopTools_ListOfShape& c);
+
+		bool get_edge_axis(const TopoDS_Edge& e, gp_Ax1& ax);
+
+		bool is_subset(const TopTools_IndexedMapOfShape& lhs, const TopTools_IndexedMapOfShape& rhs);
+
+		bool is_extrusion(const gp_Vec& v, const TopoDS_Shape& s, TopoDS_Face& base, std::pair<double, double>& interval);
+
+		int eliminate_touching_operands(double prec, const TopoDS_Shape& a, const TopTools_ListOfShape& bs, TopTools_ListOfShape& c);
+
+		TopoDS_Shape unify(const TopoDS_Shape& s, double tolerance);
+
+		bool boolean_subtraction_2d_using_builder(const TopoDS_Shape& a_input, const TopTools_ListOfShape& b_input, TopoDS_Shape& result, double eps);
+
+	}
+}
+
+#endif

src/ifcgeom_schema_agnostic/face_definition.cpp

@@ -0,0 +1,22 @@
+#include "face_definition.h"
+
+#include <TopoDS.hxx>
+#include <Geom_Line.hxx>
+#include <BRep_Tool.hxx>
+#include <TopoDS_Iterator.hxx>
+
+/* Returns whether wire conforms to a polyhedron, i.e. only edges with linear curves*/
+bool IfcGeom::util::is_polyhedron(const TopoDS_Wire & wire) {
+	double a, b;
+	TopLoc_Location l;
+
+	TopoDS_Iterator it(wire, false, false);
+	for (; it.More(); it.Next()) {
+		auto crv = BRep_Tool::Curve(TopoDS::Edge(it.Value()), l, a, b);
+		if (!crv || crv->DynamicType() != STANDARD_TYPE(Geom_Line)) {
+			return false;
+		}
+	}
+
+	return true;
+}

src/ifcgeom_schema_agnostic/face_definition.h

@@ -0,0 +1,78 @@
+/********************************************************************************
+ *                                                                              *
+ * This file is part of IfcOpenShell.                                           *
+ *                                                                              *
+ * IfcOpenShell is free software: you can redistribute it and/or modify         *
+ * it under the terms of the Lesser GNU General Public License as published by  *
+ * the Free Software Foundation, either version 3.0 of the License, or          *
+ * (at your option) any later version.                                          *
+ *                                                                              *
+ * IfcOpenShell is distributed in the hope that it will be useful,              *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of               *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the                 *
+ * Lesser GNU General Public License for more details.                          *
+ *                                                                              *
+ * You should have received a copy of the Lesser GNU General Public License     *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>.         *
+ *                                                                              *
+ ********************************************************************************/
+
+#ifndef FACE_DEFINITION_H
+#define FACE_DEFINITION_H
+
+#include <TopoDS_Wire.hxx>
+#include <Geom_Surface.hxx>
+
+#include <map>
+#include <vector>
+
+namespace IfcGeom {
+	namespace util {
+
+		/* Returns whether wire conforms to a polyhedron, i.e. only edges with linear curves*/
+		bool is_polyhedron(const TopoDS_Wire& wire);
+
+		/* A temporary structure to store the intermediate data for the face conversion */
+		class face_definition {
+		private:
+			Handle(Geom_Surface) surface_;
+			std::vector<TopoDS_Wire> wires_;
+			bool all_outer_;
+		public:
+			face_definition() : surface_(), all_outer_(false) {}
+
+			typedef std::vector<TopoDS_Wire>::const_iterator wire_it;
+
+			bool& all_outer() {
+				return all_outer_;
+			}
+
+			bool all_outer() const {
+				return all_outer_;
+			}
+
+			Handle(Geom_Surface)& surface() {
+				return surface_;
+			}
+
+			const Handle(Geom_Surface)& surface() const {
+				return surface_;
+			}
+
+			std::vector<TopoDS_Wire>& wires() {
+				return wires_;
+			}
+
+			const TopoDS_Wire& outer_wire() const {
+				return wires_.front();
+			}
+
+			std::pair<wire_it, wire_it> inner_wires() const {
+				return { wires_.begin() + 1, wires_.end() };
+			}
+		};	
+
+	}
+}
+
+#endif

src/ifcgeom_schema_agnostic/sweep_utils.cpp

@@ -0,0 +1,356 @@
+#include "sweep_utils.h"
+
+#include "../ifcparse/IfcLogger.h"
+#include "../ifcgeom_schema_agnostic/Kernel.h"
+
+#include <gp_Ax2.hxx>
+#include <Geom_Line.hxx>
+#include <Geom_Circle.hxx>
+
+#include <TopExp.hxx>
+#include <TopExp_Explorer.hxx>
+
+#include <TopoDS.hxx>
+#include <TopoDS_Vertex.hxx>
+#include <TopoDS_Compound.hxx>
+
+#include <BRep_Tool.hxx>
+#include <BRep_Builder.hxx>
+
+#include <BRepPrimAPI_MakePrism.hxx>
+#include <BRepPrimAPI_MakeRevol.hxx>
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
+#include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepOffsetAPI_MakePipeShell.hxx>
+
+bool IfcGeom::util::wire_is_c1_continuous(const TopoDS_Wire & w, double tol) {
+	// NB Note that c0 continuity is NOT checked!
+
+	TopTools_IndexedDataMapOfShapeListOfShape map;
+	TopExp::MapShapesAndAncestors(w, TopAbs_VERTEX, TopAbs_EDGE, map);
+	for (int i = 1; i <= map.Extent(); ++i) {
+		const auto& li = map.FindFromIndex(i);
+		if (li.Extent() == 2) {
+			const TopoDS_Vertex& v = TopoDS::Vertex(map.FindKey(i));
+
+			const TopoDS_Edge& e0 = TopoDS::Edge(li.First());
+			const TopoDS_Edge& e1 = TopoDS::Edge(li.Last());
+
+			double u0 = BRep_Tool::Parameter(v, e0);
+			double u1 = BRep_Tool::Parameter(v, e1);
+
+			double _, __;
+			Handle(Geom_Curve) c0 = BRep_Tool::Curve(e0, _, __);
+			Handle(Geom_Curve) c1 = BRep_Tool::Curve(e1, _, __);
+
+			gp_Pnt p;
+			gp_Vec v0, v1;
+			c0->D1(u0, p, v0);
+			c1->D1(u1, p, v1);
+
+			if (1. - std::abs(v0.Normalized().Dot(v1.Normalized())) > tol) {
+				return false;
+			}
+		}
+	}
+	return true;
+}
+
+bool IfcGeom::util::wire_to_ax(const TopoDS_Wire & wire, gp_Ax2 & directrix) {
+	gp_Pnt directrix_origin;
+	gp_Vec directrix_tangent;
+
+	TopoDS_Edge edge;
+
+	// Find first edge
+	TopoDS_Vertex v0, v1;
+	TopExp::Vertices(wire, v0, v1);
+	TopTools_IndexedDataMapOfShapeListOfShape map;
+	TopExp::MapShapesAndAncestors(wire, TopAbs_VERTEX, TopAbs_EDGE, map);
+	if (v0.IsSame(v1) && map.Contains(v0) && map.FindFromKey(v0).Extent() == 2) {
+		// Closed wire, with more than 1 edges
+		auto es = map.FindFromKey(v0);
+		auto e1 = TopoDS::Edge(es.First());
+		auto e2 = TopoDS::Edge(es.Last());
+
+		double u0, u1;
+
+		gp_Vec accum;
+
+		Handle(Geom_Curve) crv = BRep_Tool::Curve(e1, u0, u1);
+		crv->D1(TopExp::FirstVertex(e1).IsSame(v0) ? u0 : u1, directrix_origin, directrix_tangent);
+
+		accum += directrix_tangent;
+
+		crv = BRep_Tool::Curve(e2, u0, u1);
+		crv->D1(TopExp::FirstVertex(e2).IsSame(v0) ? u0 : u1, directrix_origin, directrix_tangent);
+
+		accum += directrix_tangent;
+
+		directrix_tangent = accum;
+
+	} else if (map.Contains(v0) && map.FindFromKey(v0).Extent() == 1) {
+		edge = TopoDS::Edge(map.FindFromKey(v0).First());
+
+		double u0, u1;
+		Handle(Geom_Curve) crv = BRep_Tool::Curve(edge, u0, u1);
+		crv->D1(u0, directrix_origin, directrix_tangent);
+	} else {
+		Logger::Error("Unable to locate first edge");
+		return false;
+	}
+
+	directrix = gp_Ax2(directrix_origin, directrix_tangent);
+
+	return true;
+}
+
+bool IfcGeom::util::is_single_linear_edge(const TopoDS_Wire & wire) {
+	TopExp_Explorer exp(wire, TopAbs_EDGE);
+	if (!exp.More()) {
+		return false;
+	}
+	TopoDS_Edge e = TopoDS::Edge(exp.Current());
+	exp.Next();
+	if (exp.More()) {
+		return false;
+	}
+	double u, v;
+	Handle_Geom_Curve crv = BRep_Tool::Curve(e, u, v);
+	return crv->DynamicType() == STANDARD_TYPE(Geom_Line);
+}
+
+bool IfcGeom::util::is_single_circular_edge(const TopoDS_Wire & wire) {
+	TopExp_Explorer exp(wire, TopAbs_EDGE);
+	if (!exp.More()) {
+		return false;
+	}
+	TopoDS_Edge e = TopoDS::Edge(exp.Current());
+	exp.Next();
+	if (exp.More()) {
+		return false;
+	}
+	double u, v;
+	Handle_Geom_Curve crv = BRep_Tool::Curve(e, u, v);
+	return crv->DynamicType() == STANDARD_TYPE(Geom_Circle);
+}
+
+void IfcGeom::util::process_sweep_as_extrusion(const TopoDS_Wire & wire, const TopoDS_Wire & section, TopoDS_Shape & result) {
+	TopExp_Explorer exp(wire, TopAbs_EDGE);
+	TopoDS_Edge e = TopoDS::Edge(exp.Current());
+	double u, v;
+	Handle_Geom_Curve crv = BRep_Tool::Curve(e, u, v);
+	const auto& dir = Handle(Geom_Line)::DownCast(crv)->Position().Direction();
+	// OCCT line is normalized so diff in parametric coords equals length
+	const double depth = std::abs(u - v);
+	// @todo we could be extruding the wire only when we know this is an intermediate edge.
+	TopoDS_Face face = BRepBuilderAPI_MakeFace(section).Face();
+	result = BRepPrimAPI_MakePrism(face, depth*dir).Shape();
+}
+
+void IfcGeom::util::process_sweep_as_revolution(const TopoDS_Wire & wire, const TopoDS_Wire & section, TopoDS_Shape & result) {
+	TopExp_Explorer exp(wire, TopAbs_EDGE);
+	TopoDS_Edge e = TopoDS::Edge(exp.Current());
+	double u, v;
+	Handle_Geom_Curve crv = BRep_Tool::Curve(e, u, v);
+	auto circ = Handle(Geom_Circle)::DownCast(crv);
+	// @todo we could be extruding the wire only when we know this is an intermediate edge.
+	const double depth = std::abs(u - v);
+	TopoDS_Face face = BRepBuilderAPI_MakeFace(section).Face();
+	result = BRepPrimAPI_MakeRevol(face, circ->Axis(), depth).Shape();
+}
+
+void IfcGeom::util::process_sweep_as_pipe(const TopoDS_Wire & wire, const TopoDS_Wire & section, TopoDS_Shape & result, bool force_transformed) {
+	// This tolerance is fairly high due to the linear edge substitution for small (or large radii) conical curves.
+	const bool is_continuous = wire_is_c1_continuous(wire, 1.e-2);
+	BRepOffsetAPI_MakePipeShell builder(wire);
+	builder.Add(section);
+	builder.SetTransitionMode(is_continuous || force_transformed ? BRepBuilderAPI_Transformed : BRepBuilderAPI_RightCorner);
+	try {
+		builder.Build();
+	} catch (Standard_Failure& e) {
+		// We fallback to BRepBuilderAPI_Transformed, but likely with visual artefacts.
+		if (!(is_continuous || force_transformed)) {
+			return process_sweep_as_pipe(wire, section, result, true);
+		} else {
+			throw e;
+		}
+	}
+	builder.MakeSolid();
+	result = builder.Shape();
+}
+
+void IfcGeom::util::sort_edges(const TopoDS_Wire & wire, std::vector<TopoDS_Edge>& sorted_edges) {
+	TopTools_IndexedDataMapOfShapeListOfShape map;
+	TopExp::MapShapesAndAncestors(wire, TopAbs_VERTEX, TopAbs_EDGE, map);
+
+	for (int i = 1; i <= map.Extent(); ++i) {
+		if (map.FindFromIndex(i).Extent() > 2) {
+			Logger::Warning("Self-intersecting Directrix");
+		}
+	}
+
+	std::set<TopoDS_TShape*> seen;
+
+	auto num_edges = IfcGeom::Kernel::count(wire, TopAbs_EDGE);
+
+	TopoDS_Vertex v0, v1;
+	// @todo this creates the ancestor map twice
+	TopExp::Vertices(wire, v0, v1);
+
+	bool ignore_first_equality_because_closed = v0.IsSame(v1);
+
+	// @todo this probably still does not work on a closed wire consisting of one (circular) edge.
+
+	while ((int)sorted_edges.size() < num_edges &&
+		(!v0.IsSame(v1) || ignore_first_equality_because_closed)) {
+		ignore_first_equality_because_closed = false;
+		if (!map.Contains(v0)) {
+			throw std::runtime_error("Disconnected vertex");
+		}
+		const TopTools_ListOfShape& es = map.FindFromKey(v0);
+		TopoDS_Vertex ve0, ve1;
+		TopTools_ListIteratorOfListOfShape it(es);
+		bool added = false;
+		for (; it.More(); it.Next()) {
+			const TopoDS_Edge& e = TopoDS::Edge(it.Value());
+			TopExp::Vertices(e, ve0, ve1, true);
+			if (ve0.IsSame(v0) && seen.find(&*e.TShape()) == seen.end()) {
+				sorted_edges.push_back(e);
+				v0 = ve1;
+				added = true;
+				seen.insert(&*e.TShape());
+				break;
+			}
+		}
+		if (!added) {
+			throw std::runtime_error("Disconnected edge");
+		}
+	}
+}
+
+// #939: a closed loop causes failed triangulation in 7.3 and artefacts
+// in 7.4 so we break up a closed wire into two equal parts.
+void IfcGeom::util::break_closed(const TopoDS_Wire & wire, std::vector<TopoDS_Wire>& wires) {
+	std::vector<TopoDS_Edge> sorted_edges;
+	sort_edges(wire, sorted_edges);
+
+	if (sorted_edges.size() == 1) {
+		wires.push_back(wire);
+		return;
+	}
+
+	BRep_Builder B;
+
+	wires.emplace_back();
+	B.MakeWire(wires.back());
+
+	for (size_t i = 0; i < sorted_edges.size(); ++i) {
+		if (i == sorted_edges.size() / 2) {
+			wires.emplace_back();
+			B.MakeWire(wires.back());
+		}
+
+		const auto& e = sorted_edges[i];
+		B.Add(wires.back(), e);
+	}
+}
+
+void IfcGeom::util::segment_adjacent_non_linear(const TopoDS_Wire & wire, std::vector<TopoDS_Wire>& wires) {
+	std::vector<TopoDS_Edge> sorted_edges;
+	sort_edges(wire, sorted_edges);
+
+	BRep_Builder B;
+	double u, v;
+
+	wires.emplace_back();
+	B.MakeWire(wires.back());
+
+	for (int i = 0; i < (int)sorted_edges.size() - 1; ++i) {
+		const auto& e = sorted_edges[i];
+		Handle_Geom_Curve crv = BRep_Tool::Curve(e, u, v);
+		const bool is_linear = crv->DynamicType() == STANDARD_TYPE(Geom_Line);
+
+		const auto& f = sorted_edges[i + 1];
+		crv = BRep_Tool::Curve(f, u, v);
+		const bool next_is_linear = crv->DynamicType() == STANDARD_TYPE(Geom_Line);
+
+		B.Add(wires.back(), e);
+
+		if (!is_linear && !next_is_linear) {
+			wires.emplace_back();
+			B.MakeWire(wires.back());
+		}
+	}
+
+	if (!sorted_edges.empty()) {
+		B.Add(wires.back(), sorted_edges.back());
+	}
+}
+
+// @todo make this generic for other sweeps not just swept disk
+void IfcGeom::util::process_sweep(const TopoDS_Wire & wire, double radius, TopoDS_Shape & result) {
+	std::vector<TopoDS_Wire> wires, wires_tmp;
+	segment_adjacent_non_linear(wire, wires_tmp);
+	for (auto& w : wires_tmp) {
+		break_closed(w, wires);
+	}
+
+	TopoDS_Compound C;
+	BRep_Builder B;
+	if (wires.size() > 1) {
+		B.MakeCompound(C);
+	}
+
+	for (auto& w : wires) {
+		TopoDS_Shape part;
+
+		gp_Ax2 directrix;
+		if (!wire_to_ax(w, directrix)) {
+			continue;
+		}
+		Handle(Geom_Circle) circle = new Geom_Circle(directrix, radius);
+		TopoDS_Wire section = BRepBuilderAPI_MakeWire(BRepBuilderAPI_MakeEdge(circle));
+
+		if (is_single_circular_edge(w)) {
+			process_sweep_as_revolution(w, section, part);
+		} else if (is_single_linear_edge(w)) {
+			process_sweep_as_extrusion(w, section, part);
+		} else {
+			process_sweep_as_pipe(w, section, part);
+		}
+		if (wires.size() > 1) {
+			B.Add(C, part);
+		} else {
+			result = part;
+		}
+	}
+
+	if (wires.size() > 1) {
+		result = C;
+	}
+
+	/*
+	// Eliminate Swept Surfaces?
+	result = ShapeCustom::SweptToElementary(result);
+
+	// Eliminate Trimmed Surfaces?
+	ShapeBuild_ReShape sbrs;
+	BRep_Builder b;
+	TopExp_Explorer exp(result, TopAbs_FACE);
+	for (; exp.More(); exp.Next()) {
+	const TopoDS_Face& f = TopoDS::Face(exp.Current());
+	auto S = BRep_Tool::Surface(f);
+	if (S->IsKind(STANDARD_TYPE(Geom_RectangularTrimmedSurface))) {
+	auto RTS = Handle(Geom_RectangularTrimmedSurface)::DownCast(S);
+	auto B = RTS->BasisSurface();
+	TopoDS_Shape newf = f.EmptyCopied();
+	// @todo Is it ok to assume no location?
+	b.MakeFace(TopoDS::Face(newf), B, BRep_Tool::Tolerance(f));
+	sbrs.Replace(f, newf);
+	}
+	}
+	result = sbrs.Apply(result);
+	*/
+}

src/ifcgeom_schema_agnostic/sweep_utils.h

@@ -0,0 +1,59 @@
+/********************************************************************************
+ *                                                                              *
+ * This file is part of IfcOpenShell.                                           *
+ *                                                                              *
+ * IfcOpenShell is free software: you can redistribute it and/or modify         *
+ * it under the terms of the Lesser GNU General Public License as published by  *
+ * the Free Software Foundation, either version 3.0 of the License, or          *
+ * (at your option) any later version.                                          *
+ *                                                                              *
+ * IfcOpenShell is distributed in the hope that it will be useful,              *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of               *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the                 *
+ * Lesser GNU General Public License for more details.                          *
+ *                                                                              *
+ * You should have received a copy of the Lesser GNU General Public License     *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>.         *
+ *                                                                              *
+ ********************************************************************************/
+
+#ifndef SWEEP_UTILS_H
+#define SWEEP_UTILS_H
+
+#include <TopoDS_Wire.hxx>
+#include <TopoDS_Edge.hxx>
+
+#include <vector>
+
+namespace IfcGeom {
+	namespace util {
+
+		bool wire_is_c1_continuous(const TopoDS_Wire& w, double tol);
+
+		bool wire_to_ax(const TopoDS_Wire& wire, gp_Ax2& directrix);
+
+		bool is_single_linear_edge(const TopoDS_Wire& wire);
+
+		bool is_single_circular_edge(const TopoDS_Wire& wire);
+
+		void process_sweep_as_extrusion(const TopoDS_Wire& wire, const TopoDS_Wire& section, TopoDS_Shape& result);
+
+		void process_sweep_as_revolution(const TopoDS_Wire& wire, const TopoDS_Wire& section, TopoDS_Shape& result);
+
+		void process_sweep_as_pipe(const TopoDS_Wire& wire, const TopoDS_Wire& section, TopoDS_Shape& result, bool force_transformed = false);
+
+		void sort_edges(const TopoDS_Wire& wire, std::vector<TopoDS_Edge>& sorted_edges);
+
+
+		// #939: a closed loop causes failed triangulation in 7.3 and artefacts
+		// in 7.4 so we break up a closed wire into two equal parts.
+		void break_closed(const TopoDS_Wire& wire, std::vector<TopoDS_Wire>& wires);
+
+		void segment_adjacent_non_linear(const TopoDS_Wire& wire, std::vector<TopoDS_Wire>& wires);
+
+		// @todo make this generic for other sweeps not just swept disk
+		void process_sweep(const TopoDS_Wire& wire, double radius, TopoDS_Shape& result);
+	}
+}
+
+#endif

src/ifcgeom_schema_agnostic/wire_builder.cpp

@@ -0,0 +1,226 @@
+#include "wire_builder.h"
+
+#include "../ifcparse/IfcLogger.h"
+#include "../ifcgeom_schema_agnostic/Kernel.h"
+
+#include <TopExp.hxx>
+#include <TopoDS.hxx>
+#include <BRep_Tool.hxx>
+#include <BRep_Builder.hxx>
+#include <ShapeBuild_ReShape.hxx>
+#include <GC_MakeCircle.hxx>
+#include <BRepBuilderAPI_MakeEdge.hxx>
+#include <Geom_Line.hxx>
+#include <Geom_Circle.hxx>
+#include <GeomAdaptor_Curve.hxx>
+
+// Returns the first edge of a wire
+TopoDS_Edge IfcGeom::util::first_edge(const TopoDS_Wire & w) {
+	TopoDS_Vertex v1, v2;
+	TopExp::Vertices(w, v1, v2);
+	TopTools_IndexedDataMapOfShapeListOfShape wm;
+	TopExp::MapShapesAndAncestors(w, TopAbs_VERTEX, TopAbs_EDGE, wm);
+	return TopoDS::Edge(wm.FindFromKey(v1).First());
+}
+
+// Returns new wire with the edge replaced by a linear edge with the vertex v moved to p
+TopoDS_Wire IfcGeom::util::adjust(const TopoDS_Wire & w, const TopoDS_Vertex & v, const gp_Pnt & p) {
+	TopTools_IndexedDataMapOfShapeListOfShape map;
+	TopExp::MapShapesAndAncestors(w, TopAbs_VERTEX, TopAbs_EDGE, map);
+
+	bool all_linear = true, single_circle = false, first = true;
+
+	const TopTools_ListOfShape& edges = map.FindFromKey(v);
+	TopTools_ListIteratorOfListOfShape it(edges);
+	for (; it.More(); it.Next()) {
+		const TopoDS_Edge& e = TopoDS::Edge(it.Value());
+		double _, __;
+		Handle(Geom_Curve) crv = BRep_Tool::Curve(e, _, __);
+		const bool is_line = crv->DynamicType() == STANDARD_TYPE(Geom_Line);
+		const bool is_circle = crv->DynamicType() == STANDARD_TYPE(Geom_Circle);
+		all_linear = all_linear && is_line;
+		single_circle = first && is_circle;
+	}
+
+	if (all_linear) {
+		BRep_Builder b;
+		TopoDS_Vertex v2;
+		b.MakeVertex(v2, p, BRep_Tool::Tolerance(v));
+
+		ShapeBuild_ReShape reshape;
+		reshape.Replace(v.Oriented(TopAbs_FORWARD), v2);
+
+		return TopoDS::Wire(reshape.Apply(w));
+	} else if (single_circle) {
+		TopoDS_Vertex v1, v2;
+		TopExp::Vertices(w, v1, v2);
+
+		gp_Pnt p1, p2, p3;
+		p1 = v.IsEqual(v1) ? p : BRep_Tool::Pnt(v1);
+		p3 = v.IsEqual(v2) ? p : BRep_Tool::Pnt(v2);
+
+		double a, b;
+		Handle(Geom_Curve) crv = BRep_Tool::Curve(TopoDS::Edge(edges.First()), a, b);
+		crv->D0((a + b) / 2., p2);
+
+		GC_MakeCircle mc(p1, p2, p3);
+		if (!mc.IsDone()) {
+			throw IfcGeom::geometry_exception("Failed to adjust circle");
+		}
+
+		TopoDS_Edge edge = BRepBuilderAPI_MakeEdge(mc.Value(), p1, p3).Edge();
+		BRepBuilderAPI_MakeWire builder;
+		builder.Add(edge);
+		return builder.Wire();
+	} else {
+		throw IfcGeom::geometry_exception("Unexpected wire to adjust");
+	}
+}
+
+double IfcGeom::util::deflection_for_approximating_circle(double radius, double param) {
+	return -radius * std::cos(1. / 2. * param) * std::cos(param) - radius * std::sin(1. / 2. * param) * std::sin(param) + radius;
+}
+
+bool IfcGeom::util::create_edge_over_curve_with_log_messages(const Handle_Geom_Curve & crv, const double eps, const gp_Pnt & p1, const gp_Pnt & p2, TopoDS_Edge & result) {
+	if (crv->IsClosed() && p1.Distance(p2) <= eps) {
+		BRepBuilderAPI_MakeEdge me(crv);
+		if (me.IsDone()) {
+			result = me.Edge();
+			return true;
+		} else {
+			return false;
+		}
+	}
+
+	BRep_Builder builder;
+	TopoDS_Vertex v1, v2;
+	/// @todo project first and emit warnings accordingly
+	builder.MakeVertex(v1, p1, eps);
+	builder.MakeVertex(v2, p2, eps);
+
+	BRepBuilderAPI_MakeEdge me(crv, v1, v2);
+	if (!me.IsDone()) {
+		const double eps2 = eps * eps;
+		if (me.Error() == BRepBuilderAPI_PointProjectionFailed) {
+			GeomAdaptor_Curve GAC(crv);
+			const gp_Pnt* ps[2] = { &p1, &p2 };
+			for (int i = 0; i < 2; ++i) {
+				Extrema_ExtPC extrema(*ps[i], GAC);
+				if (extrema.IsDone()) {
+					int n = extrema.NbExt();
+					double dmin = std::numeric_limits<double>::infinity();
+					for (int j = 1; j <= n; j++) {
+						const double d = extrema.SquareDistance(j);
+						if (d < dmin) {
+							dmin = d;
+						}
+					}
+					if (dmin == std::numeric_limits<double>::infinity()) {
+						Logger::Error("No extrema for point");
+					} else if (dmin > eps2) {
+						Logger::Error("Distance of " + boost::lexical_cast<std::string>(std::sqrt(dmin)) + " exceeds tolerance");
+					}
+				} else {
+					Logger::Error("Failed to calculate extrema for point");
+				}
+			}
+		}
+		return false;
+	}
+	result = me.Edge();
+	return true;
+}
+
+void IfcGeom::util::wire_builder::operator()(const TopoDS_Shape& a) {
+	const TopoDS_Wire& w = TopoDS::Wire(a);
+	if (override_next_) {
+		override_next_ = false;
+		TopoDS_Edge e = first_edge(w);
+		mw_.Add(adjust(w, TopExp::FirstVertex(e, true), next_override_));
+	} else {
+		mw_.Add(w);
+	}
+}
+
+void IfcGeom::util::wire_builder::operator()(const TopoDS_Shape& a, const TopoDS_Shape& b, bool last) {
+	TopoDS_Wire w1 = TopoDS::Wire(a);
+	const TopoDS_Wire& w2 = TopoDS::Wire(b);
+
+	if (override_next_) {
+		override_next_ = false;
+		TopoDS_Edge e = first_edge(w1);
+		w1 = adjust(w1, TopExp::FirstVertex(e, true), next_override_);
+	}
+
+	TopoDS_Vertex w11, w12, w21, w22;
+	TopExp::Vertices(w1, w11, w12);
+	TopExp::Vertices(w2, w21, w22);
+
+	gp_Pnt p1 = BRep_Tool::Pnt(w12);
+	gp_Pnt p2 = BRep_Tool::Pnt(w21);
+
+	double dist = p1.Distance(p2);
+
+	// Distance is within tolerance, this is fine
+	if (dist < p_) {
+		mw_.Add(w1);
+		goto check;
+	}
+
+	// Distance is too large for attempting to move end points, add intermediate edge
+	if (dist > 1000. * p_) {
+		mw_.Add(w1);
+		mw_.Add(BRepBuilderAPI_MakeEdge(p1, p2));
+		Logger::Warning("Added additional segment to close gap with length " + boost::lexical_cast<std::string>(dist) + " to:", inst_);
+		goto check;
+	}
+
+	{
+		TopTools_IndexedDataMapOfShapeListOfShape wmap1, wmap2;
+
+		// Find edges connected to end- and begin vertex
+		TopExp::MapShapesAndAncestors(w1, TopAbs_VERTEX, TopAbs_EDGE, wmap1);
+		TopExp::MapShapesAndAncestors(w2, TopAbs_VERTEX, TopAbs_EDGE, wmap2);
+
+		const TopTools_ListOfShape& last_edges = wmap1.FindFromKey(w12);
+		const TopTools_ListOfShape& first_edges = wmap2.FindFromKey(w21);
+
+		double _, __;
+		if (last_edges.Extent() == 1 && first_edges.Extent() == 1) {
+			Handle(Geom_Curve) c1 = BRep_Tool::Curve(TopoDS::Edge(last_edges.First()), _, __);
+			Handle(Geom_Curve) c2 = BRep_Tool::Curve(TopoDS::Edge(first_edges.First()), _, __);
+
+			const bool is_line1 = c1->DynamicType() == STANDARD_TYPE(Geom_Line);
+			const bool is_line2 = c2->DynamicType() == STANDARD_TYPE(Geom_Line);
+
+			const bool is_circle1 = c1->DynamicType() == STANDARD_TYPE(Geom_Circle);
+			const bool is_circle2 = c2->DynamicType() == STANDARD_TYPE(Geom_Circle);
+
+			// Preferably adjust the segment that is linear
+			if (is_line1 || (is_circle1 && !is_line2)) {
+				mw_.Add(adjust(w1, w12, p2));
+				Logger::Notice("Adjusted edge end-point with distance " + boost::lexical_cast<std::string>(dist) + " on:", inst_);
+			} else if ((is_line2 || is_circle2) && !last) {
+				mw_.Add(w1);
+				override_next_ = true;
+				next_override_ = p1;
+				Logger::Notice("Adjusted edge end-point with distance " + boost::lexical_cast<std::string>(dist) + " on:", inst_);
+			} else {
+				// In all other cases an edge is added
+				mw_.Add(w1);
+				mw_.Add(BRepBuilderAPI_MakeEdge(p1, p2));
+				Logger::Warning("Added additional segment to close gap with length " + boost::lexical_cast<std::string>(dist) + " to:", inst_);
+			}
+		} else {
+			Logger::Error("Internal error, inconsistent wire segments", inst_);
+			mw_.Add(w1);
+		}
+	}
+
+check:
+	if (mw_.Error() == BRepBuilderAPI_NonManifoldWire) {
+		Logger::Error("Non-manifold curve segments:", inst_);
+	} else if (mw_.Error() == BRepBuilderAPI_DisconnectedWire) {
+		Logger::Error("Failed to join curve segments:", inst_);
+	}
+}

src/ifcgeom_schema_agnostic/wire_builder.h

@@ -0,0 +1,105 @@
+/********************************************************************************
+ *                                                                              *
+ * This file is part of IfcOpenShell.                                           *
+ *                                                                              *
+ * IfcOpenShell is free software: you can redistribute it and/or modify         *
+ * it under the terms of the Lesser GNU General Public License as published by  *
+ * the Free Software Foundation, either version 3.0 of the License, or          *
+ * (at your option) any later version.                                          *
+ *                                                                              *
+ * IfcOpenShell is distributed in the hope that it will be useful,              *
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of               *
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the                 *
+ * Lesser GNU General Public License for more details.                          *
+ *                                                                              *
+ * You should have received a copy of the Lesser GNU General Public License     *
+ * along with this program. If not, see <http://www.gnu.org/licenses/>.         *
+ *                                                                              *
+ ********************************************************************************/
+
+#ifndef WIRE_BUILDER_H
+#define WIRE_BUILDER_H
+
+#include "../ifcparse/IfcBaseClass.h"
+
+#include <Geom_Curve.hxx>
+
+#include <TopoDS_Vertex.hxx>
+#include <TopoDS_Edge.hxx>
+#include <TopoDS_Wire.hxx>
+
+#include <Extrema_ExtPC.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
+
+namespace IfcGeom {
+	namespace util {
+		// Returns the first edge of a wire
+		TopoDS_Edge first_edge(const TopoDS_Wire& w);
+
+		// Returns new wire with the edge replaced by a linear edge with the vertex v moved to p
+		TopoDS_Wire adjust(const TopoDS_Wire& w, const TopoDS_Vertex& v, const gp_Pnt& p);
+
+		// A wrapper around BRepBuilderAPI_MakeWire that makes sure segments are connected either by moving end points or by adding intermediate segments
+		class wire_builder {
+		private:
+			BRepBuilderAPI_MakeWire mw_;
+			double p_;
+			bool override_next_;
+			gp_Pnt next_override_;
+			const IfcUtil::IfcBaseClass* inst_;
+
+		public:
+			wire_builder(double p, const IfcUtil::IfcBaseClass* inst = 0) : p_(p), override_next_(false), inst_(inst) {}
+
+			void operator()(const TopoDS_Shape& a);
+
+			void operator()(const TopoDS_Shape& a, const TopoDS_Shape& b, bool last);
+
+			const TopoDS_Wire& wire() { return mw_.Wire(); }
+		};
+
+		template <typename Fn>
+		void shape_pair_enumerate(TopTools_ListIteratorOfListOfShape& it, Fn& fn, bool closed) {
+			bool is_first = true;
+			TopoDS_Shape first, previous, current;
+			for (; it.More(); it.Next(), is_first = false) {
+				current = it.Value();
+				if (is_first) {
+					first = current;
+				} else {
+					fn(previous, current, false);
+				}
+				previous = current;
+			}
+			if (closed) {
+				fn(current, first, true);
+			} else {
+				fn(current);
+			}
+		}
+
+		/*
+		Below is code to deduce the formula below in SageMath
+
+		| R, b = var('R b')
+		|
+		| Bxy = R * cos(b), R * sin(b)
+		| Cxy = R * cos(b/2), R * sin(b/2)
+		|
+		| def dot(v, w):
+		|     return v[0] * w[0] + v[1] * w[1]
+		|
+		| def norm(v):
+		|     l = sqrt(v[0]^2 + v[1]^2)
+		|     return v[0] / l, v[1] / l
+		|
+		| (R - R*dot(norm(Cxy), norm(Bxy))).full_simplify()
+		*/
+
+		double deflection_for_approximating_circle(double radius, double param);
+
+		bool create_edge_over_curve_with_log_messages(const Handle_Geom_Curve& crv, const double eps, const gp_Pnt& p1, const gp_Pnt& p2, TopoDS_Edge& result);
+	}
+}
+
+#endif