IfcOpenShell/src/ifcgeom/kernels/opencascade/loop.cpp

#include "OpenCascadeKernel.h"
#include "wire_builder.h"

#include <Geom_Line.hxx>
#include <Geom_Circle.hxx>
#include <Geom_Ellipse.hxx>
#include <BRepAdaptor_CompCurve.hxx>
#include <Approx_Curve3d.hxx>
#include <ShapeFix_ShapeTolerance.hxx>
#include <BRepBuilderAPI_MakeWire.hxx>
#include <gp_Pnt.hxx>
#include <TColgp_Array1OfPnt.hxx>
#include <TColStd_Array1OfReal.hxx>
#include <TColStd_Array1OfInteger.hxx>
#include <Geom_BSplineCurve.hxx>

#include <Standard_Version.hxx>
#if OCC_VERSION_HEX < 0x70600
#include <BRepAdaptor_HCompCurve.hxx>
#endif

using namespace ifcopenshell::geometry;
using namespace ifcopenshell::geometry::kernels;
using namespace IfcGeom;
using namespace IfcGeom::util;

namespace {
	typedef boost::variant<Handle(Geom_Curve), TopoDS_Wire> curve_creation_visitor_result_type;
	curve_creation_visitor_result_type convert_curve(OpenCascadeKernel* kernel, const taxonomy::item* curve);

	struct curve_creation_visitor {
		OpenCascadeKernel* kernel;
		curve_creation_visitor_result_type result;

		curve_creation_visitor_result_type operator()(const taxonomy::bspline_curve& bc) {

			const bool is_rational = !!bc.weights;

			TColgp_Array1OfPnt      Poles(0, bc.control_points.size() - 1);
			TColStd_Array1OfReal    Weights(0, bc.control_points.size() - 1);
			TColStd_Array1OfReal    Knots(0, (int)bc.knots.size() - 1);
			TColStd_Array1OfInteger Mults(0, (int)bc.knots.size() - 1);
			Standard_Integer        Degree = bc.degree;
			Standard_Boolean        Periodic = false;
			// @tfk: it appears to be wrong to expect a period curve when the curve is closed, see #586
			// Standard_Boolean     Periodic = l->ClosedCurve();

			int i;

			if (is_rational) {
				i = 0;
				for (auto it = bc.weights->begin(); it != bc.weights->end(); ++it, ++i) {
					Weights(i) = *it;
				}
			}

			i = 0;
			for (auto it = bc.control_points.begin(); it != bc.control_points.end(); ++it, ++i) {
				Poles(i) = OpenCascadeKernel::convert_xyz<gp_Pnt>(*it);
			}

			i = 0;
			for (auto it = bc.multiplicities.begin(); it != bc.multiplicities.end(); ++it, ++i) {
				Mults(i) = *it;
			}

			i = 0;
			for (auto it = bc.knots.begin(); it != bc.knots.end(); ++it, ++i) {
				Knots(i) = *it;
			}

			if (is_rational) {
				return result = Handle(Geom_Curve)(new Geom_BSplineCurve(Poles, Weights, Knots, Mults, Degree, Periodic));
			} else {
				return result = Handle(Geom_Curve)(new Geom_BSplineCurve(Poles, Knots, Mults, Degree, Periodic));
			}
		}

		curve_creation_visitor_result_type operator()(const taxonomy::line& l) {
			const auto& m = l.matrix.ccomponents();
			return result = Handle(Geom_Curve)(new Geom_Line(OpenCascadeKernel::convert_xyz2<gp_Pnt>(m.col(3)), OpenCascadeKernel::convert_xyz2<gp_Dir>(m.col(0))));
		}

		curve_creation_visitor_result_type operator()(const taxonomy::circle& c) {
			const auto& m = c.matrix.ccomponents();
			return result = Handle(Geom_Curve)(new Geom_Circle(gp_Ax2(OpenCascadeKernel::convert_xyz2<gp_Pnt>(m.col(3)), OpenCascadeKernel::convert_xyz2<gp_Dir>(m.col(2)), OpenCascadeKernel::convert_xyz2<gp_Dir>(m.col(0))), c.radius));
		}

		curve_creation_visitor_result_type operator()(const taxonomy::ellipse& e) {
			const auto& m = e.matrix.ccomponents();
			return result = Handle(Geom_Curve)(new Geom_Ellipse(gp_Ax2(OpenCascadeKernel::convert_xyz2<gp_Pnt>(m.col(3)), OpenCascadeKernel::convert_xyz2<gp_Dir>(m.col(2)), OpenCascadeKernel::convert_xyz2<gp_Dir>(m.col(0))), e.radius, e.radius2));
		}

		curve_creation_visitor_result_type operator()(const taxonomy::loop& l) {
			TopoDS_Wire wire;
			kernel->convert(&l, wire);
			return result = wire;
		}

		curve_creation_visitor_result_type operator()(const taxonomy::edge& e) {
			// @todo for polyloops/-lines we should probably construct edges based on correct oriented TopoDS_Vertex instead.

			if (e.start.which() != e.end.which()) {
				throw std::runtime_error("Different trim types not supported");
			}

			TopoDS_Edge E;
			if (e.basis) {
				auto crv_or_wire = convert_curve(kernel, e.basis);
				Handle(Geom_Curve) curve;
				if (crv_or_wire.which() == 0) {
					curve = boost::get<Handle(Geom_Curve)>(crv_or_wire);
				} else {
					// @todo
					const double precision_ = 1.e-5;
					Logger::Warning("Approximating BasisCurve due to possible discontinuities", e.instance);
					const auto& w = boost::get<TopoDS_Wire>(crv_or_wire);
#if OCC_VERSION_HEX < 0x70600
					BRepAdaptor_CompCurve cc(w, true);
					Handle(Adaptor3d_HCurve) hcc = Handle(Adaptor3d_HCurve)(new BRepAdaptor_HCompCurve(cc));
#else
					auto hcc = new BRepAdaptor_CompCurve(w, true);
#endif
					// @todo, arbitrary numbers here, note they cannot be too high as contiguous memory is allocated based on them.
					Approx_Curve3d approx(hcc, precision_, GeomAbs_C0, 10, 10);
					curve = approx.Curve();
				}

				const bool reversed = !((taxonomy::geom_item*)e.basis)->orientation.get_value_or(true);
				const bool is_conic = e.basis->kind() == taxonomy::ELLIPSE || e.basis->kind() == taxonomy::CIRCLE;

				// @todo, copy over logic from previous IfcTrimmedCurve handling
				if (e.start.which() == 0) {
					auto p1 = OpenCascadeKernel::convert_xyz<gp_Pnt>(boost::get<taxonomy::point3>(e.start));
					auto p2 = OpenCascadeKernel::convert_xyz<gp_Pnt>(boost::get<taxonomy::point3>(e.end));

					if (reversed) {
						std::swap(p1, p2);
					}

					E = BRepBuilderAPI_MakeEdge(curve, p1, p2).Edge();
				} else {
					auto v1 = boost::get<double>(e.start);
					auto v2 = boost::get<double>(e.end);

					if (reversed) {
						std::swap(v1, v2);
					}

					if (is_conic && ALMOST_THE_SAME(fmod(v2 - v1, M_PI*2.), 0.)) {
						E = BRepBuilderAPI_MakeEdge(curve).Edge();
					} else {
						E = BRepBuilderAPI_MakeEdge(curve, v1, v2).Edge();
					}
				}

				if (reversed) {
					E.Reverse();
				}
			} else {
				if (e.start.which() != 0) {
					throw std::runtime_error("Non-cartesian trim on edge without curve");
				}
				auto p1 = OpenCascadeKernel::convert_xyz<gp_Pnt>(boost::get<taxonomy::point3>(e.start));
				auto p2 = OpenCascadeKernel::convert_xyz<gp_Pnt>(boost::get<taxonomy::point3>(e.end));

				E = BRepBuilderAPI_MakeEdge(p1, p2).Edge();
			}
			BRep_Builder B;
			TopoDS_Wire W;
			B.MakeWire(W);
			B.Add(W, E);
			return result = W;
		}

		curve_creation_visitor_result_type operator()(const taxonomy::offset_curve& l) {
			// @todo
			throw std::runtime_error("Offset curves not supported as part of loop");
		}
	};

	curve_creation_visitor_result_type convert_curve(OpenCascadeKernel* kernel, const taxonomy::item* curve) {
		curve_creation_visitor v{ kernel };
		if (dispatch_curve_creation<curve_creation_visitor, 0>::dispatch(curve, v)) {
			return v.result;
		} else {
			throw std::runtime_error("No curve created");
		}
	}
}

bool OpenCascadeKernel::convert(const taxonomy::loop* loop, TopoDS_Wire& wire) {
	auto segments = loop->children_as<taxonomy::edge>();

	TopTools_ListOfShape converted_segments;

	for (auto& segment : segments) {
		auto segment_wire = boost::get<TopoDS_Wire>(convert_curve(this, segment));

#ifdef IFOPSH_DEBUG
		std::ostringstream o;
		segment->print(o);
		TopoDS_Vertex v0, v1;
		TopExp::Vertices(segment_wire, v0, v1);
		gp_Pnt p0 = BRep_Tool::Pnt(v0);
		gp_Pnt p1 = BRep_Tool::Pnt(v1);
		o << "p0 " << p0.X() << " " << p0.Y() << " " << p0.Z() << std::endl;
		o << "p1 " << p1.X() << " " << p1.Y() << " " << p1.Z() << std::endl;
		auto o_str = o.str();
		std::wcout << o_str.c_str() << std::endl;
#endif

		if (!segment->orientation_2.get_value_or(true)) {
			segment_wire.Reverse();
		}

		ShapeFix_ShapeTolerance FTol;
		FTol.SetTolerance(segment_wire, precision_, TopAbs_WIRE);

		converted_segments.Append(segment_wire);
	}

	if (converted_segments.Extent() == 0) {
		Logger::Message(Logger::LOG_ERROR, "No segment succesfully converted:", loop->instance);
		return false;
	}

	BRepBuilderAPI_MakeWire w;
	TopoDS_Vertex wire_first_vertex, wire_last_vertex, edge_first_vertex, edge_last_vertex;

	TopTools_ListIteratorOfListOfShape it(converted_segments);

	/*
	@todo
	IfcEntityList::ptr profile = l->data().getInverse(&IfcSchema::IfcProfileDef::Class(), -1);
	const bool force_close = profile && profile->size() > 0;
	*/
	const bool force_close = false;

	wire_builder bld(precision_, loop->instance ? loop->instance->as<IfcUtil::IfcBaseEntity>() : nullptr);
	shape_pair_enumerate(it, bld, force_close);
	wire = bld.wire();

	return true;
}