diff --git a/src/ifcconvert/IfcConvert.cpp b/src/ifcconvert/IfcConvert.cpp index 845ad5d2c1..060b75e4c8 100644 --- a/src/ifcconvert/IfcConvert.cpp +++ b/src/ifcconvert/IfcConvert.cpp @@ -1182,6 +1182,13 @@ namespace latebound_access { #include #include +#include +#include +#include +#include +#include +namespace SMS = CGAL::Surface_mesh_simplification; + template T enlarge(const T& t, double d = 1.e-5) { T::NT min[3]; @@ -1191,6 +1198,7 @@ T enlarge(const T& t, double d = 1.e-5) { max[i] = t.max_coord(i) + d; } return T(min, max, t.handle()); + // return T(min, max); } int convert_to_nef(cgal_shape_t& shape, CGAL::Nef_polyhedron_3& result) { @@ -1227,9 +1235,576 @@ int convert_to_nef(cgal_shape_t& shape, CGAL::Nef_polyhedron_3& result) return 0; } +namespace { + // Can be used to convert polyhedron from exact to inexact and vice-versa + template + struct Copy_polyhedron_to + : public CGAL::Modifier_base { + Copy_polyhedron_to(const Polyhedron_input& in_poly) + : in_poly(in_poly) {} + + void operator()(typename Polyhedron_output::HalfedgeDS& out_hds) { + typedef typename Polyhedron_output::HalfedgeDS Output_HDS; + typedef typename Polyhedron_input::HalfedgeDS Input_HDS; + + CGAL::Polyhedron_incremental_builder_3 builder(out_hds); + + typedef typename Polyhedron_input::Vertex_const_iterator Vertex_const_iterator; + typedef typename Polyhedron_input::Facet_const_iterator Facet_const_iterator; + typedef typename Polyhedron_input::Halfedge_around_facet_const_circulator HFCC; + + builder.begin_surface(in_poly.size_of_vertices(), + in_poly.size_of_facets(), + in_poly.size_of_halfedges()); + + for (Vertex_const_iterator + vi = in_poly.vertices_begin(), end = in_poly.vertices_end(); + vi != end; ++vi) { + typename Polyhedron_output::Point_3 p(::CGAL::to_double(vi->point().x()), + ::CGAL::to_double(vi->point().y()), + ::CGAL::to_double(vi->point().z())); + builder.add_vertex(p); + } + + typedef CGAL::Inverse_index Index; + Index index(in_poly.vertices_begin(), in_poly.vertices_end()); + + for (Facet_const_iterator + fi = in_poly.facets_begin(), end = in_poly.facets_end(); + fi != end; ++fi) { + HFCC hc = fi->facet_begin(); + HFCC hc_end = hc; + builder.begin_facet(); + do { + builder.add_vertex_to_facet(index[hc->vertex()]); + ++hc; + } while (hc != hc_end); + builder.end_facet(); + } + builder.end_surface(); + } // end operator()(..) + private: + const Polyhedron_input& in_poly; + }; // end Copy_polyhedron_to<> + + template + void poly_copy(Poly_B& poly_b, const Poly_A& poly_a) { + poly_b.clear(); + Copy_polyhedron_to modifier(poly_a); + poly_b.delegate(modifier); + } + +} + +namespace { + // The following is a Visitor that keeps track of the simplification process. +// In this example the progress is printed real-time and a few statistics are +// recorded (and printed in the end). +// + struct Stats { + Stats() + : collected(0) + , processed(0) + , collapsed(0) + , non_collapsable(0) + , cost_uncomputable(0) + , placement_uncomputable(0) {} + + std::size_t collected; + std::size_t processed; + std::size_t collapsed; + std::size_t non_collapsable; + std::size_t cost_uncomputable; + std::size_t placement_uncomputable; + }; + struct My_visitor : SMS::Edge_collapse_visitor_base>> { + My_visitor(Stats* s) : stats(s) {} + // Called during the collecting phase for each edge collected. + void OnCollected(Profile const&, boost::optional const&) { + ++stats->collected; + std::wcerr << "\rEdges collected: " << stats->collected << std::flush; + } + + // Called during the processing phase for each edge selected. + // If cost is absent the edge won't be collapsed. + void OnSelected(Profile const& + , boost::optional cost + , std::size_t initial + , std::size_t current + ) { + ++stats->processed; + if (!cost) + ++stats->cost_uncomputable; + + if (current == initial) + std::wcerr << "\n" << std::flush; + std::wcerr << "\r" << current << std::flush; + } + + // Called during the processing phase for each edge being collapsed. + // If placement is absent the edge is left uncollapsed. + void OnCollapsing(Profile const& + , boost::optional placement + ) { + if (!placement) + ++stats->placement_uncomputable; + } + + // Called for each edge which failed the so called link-condition, + // that is, which cannot be collapsed because doing so would + // turn the surface mesh into a non-manifold. + void OnNonCollapsable(Profile const&) { + ++stats->non_collapsable; + } + + // Called after each edge has been collapsed + void OnCollapsed(Profile const&, vertex_descriptor) { + ++stats->collapsed; + } + + Stats* stats; + }; + +} + +namespace { + template + T approx_normalized(const T& t) { + return t * (1. / Kernel_::FT(CGAL::sqrt(CGAL::to_double(t.squared_length())))); + } +} + +#include +#include +#include +#include + +namespace { + template + struct Build_Offset : public CGAL::Modifier_base { + std::list input; + + void operator()(HDS& hds) { + // Postcondition: hds is a valid polyhedral surface. + CGAL::Polyhedron_incremental_builder_3 B(hds); + + int Nv = 0, Nf = 0; + for (auto& f : input) { + Nv += 3; + Nf += 1; + } + + B.begin_surface(Nv, Nf); + + for (auto& f : input) { + auto p0 = f->facet_begin()->vertex()->point(); + auto p1 = f->facet_begin()->next()->vertex()->point(); + auto p2 = f->facet_begin()->next()->next()->vertex()->point(); + + auto O = CGAL::centroid(p0, p1, p2); + + Kernel_::Point_3* p012[3] = { &p0, &p1, &p2 }; + for (int i = 0; i < 3; ++i) { + *p012[i] = CGAL::ORIGIN + (((*(p012[i])) - CGAL::ORIGIN) + ((*(p012[i])) - O)); + B.add_vertex(*p012[i]); + } + } + + Nv = 0; + for (int i = 0; i < Nf; ++i) { + B.begin_facet(); + B.add_vertex_to_facet(Nv++); + B.add_vertex_to_facet(Nv++); + B.add_vertex_to_facet(Nv++); + B.end_facet(); + } + + B.end_surface(); + } + }; + + template + std::list connected_faces(cgal_shape_t::Facet_handle& f, const Ts& excluded) { + std::set fs = { f }; + + std::function process; + process = [&fs, &process, &excluded](cgal_shape_t::Facet_handle& f) { + cgal_shape_t::Halfedge_around_facet_circulator circ = f->facet_begin(), end(circ); + do { + auto ff = circ->opposite()->facet(); + if (excluded.find(ff) == excluded.end()) { + auto p = fs.insert(ff); + if (p.second) { + process(ff); + } + } + } while (++circ != end); + }; + + process(f); + return std::list(fs.begin(), fs.end()); + } + + template + struct Builder_With_Map : public CGAL::Modifier_base { + std::list input; + std::map mapping; + + void operator()(HDS& hds) { + // Postcondition: hds is a valid polyhedral surface. + CGAL::Polyhedron_incremental_builder_3 B(hds); + + std::set used_points; + + for (auto& f : input) { + cgal_shape_t::Halfedge_around_facet_circulator circ = f->facet_begin(), end(circ); + do { + auto P = circ->vertex()->point(); + auto it = mapping.find(P); + if (it == mapping.end()) { + std::wcout << "WARNING unprojected point :(" << std::endl; + } else { + P = it->second; + } + used_points.insert(P); + } while (++circ != end); + } + + B.begin_surface(used_points.size(), input.size()); + + for (auto& p : used_points) { + B.add_vertex(p); + } + + for (auto& f : input) { + B.begin_facet(); + cgal_shape_t::Halfedge_around_facet_circulator circ = f->facet_begin(), end(circ); + do { + auto it = used_points.find(circ->vertex()->point()); + B.add_vertex_to_facet(std::distance(used_points.begin(), it)); + } while (++circ != end); + + B.end_facet(); + } + + B.end_surface(); + } + }; +} + +namespace { + template + T edge_collapse(T polyhedron) { + typedef CGAL::Simple_cartesian simple; + CGAL::Polyhedron_3 simple_poly; + poly_copy(simple_poly, polyhedron); + + // flattening from a thin box to a plane is not valid in edge_collapse() + Stats stats; + My_visitor vis(&stats); + SMS::Edge_length_cost elc; + SMS::Edge_length_stop_predicate stop(1.e-3); + + int r = SMS::edge_collapse(simple_poly, stop, + CGAL::parameters::vertex_index_map(get(CGAL::vertex_external_index, simple_poly)) + .halfedge_index_map(get(CGAL::halfedge_external_index, simple_poly)) + .visitor(vis) + .get_cost(elc) + ); + + std::wcout << "Removed: " << r << std::endl; + + T result; + poly_copy(result, simple_poly); + return result; + } + + + double facet_area(const cgal_shape_t::Facet_handle& f) { + auto p0 = f->facet_begin()->vertex()->point(); + auto p1 = f->facet_begin()->next()->vertex()->point(); + auto p2 = f->facet_begin()->next()->next()->vertex()->point(); + return std::sqrt(CGAL::to_double(CGAL::cross_product(p0 - p1, p2 - p1).squared_length())); + } + + void dump_facet(const cgal_shape_t::Facet_handle& f) { + auto p0 = f->facet_begin()->vertex()->point(); + auto p1 = f->facet_begin()->next()->vertex()->point(); + auto p2 = f->facet_begin()->next()->next()->vertex()->point(); + auto V = CGAL::cross_product(p0 - p1, p2 - p1); + auto d = std::sqrt(CGAL::to_double(V.squared_length())); + if (d > 1.e-20) { + V /= d; + } + + std::ostringstream oss; + oss.precision(8); + oss << "Facet with area " << facet_area(f) << " and normal (" + << CGAL::to_double(V.cartesian(0)) << " " << CGAL::to_double(V.cartesian(1)) << " " + << CGAL::to_double(V.cartesian(2)) << ")"; + + auto osss = oss.str(); + std::wcout << osss.c_str() << std::endl; + } + + struct remove_thickness { + typedef Kernel_::Point_3 Point; + typedef Kernel_::Plane_3 Plane; + typedef Kernel_::Vector_3 Vector; + typedef Kernel_::Segment_3 Segment; + typedef Kernel_::Ray_3 Ray; + + typedef CGAL::Polyhedron_3 Polyhedron; + typedef CGAL::AABB_face_graph_triangle_primitive Primitive; + typedef CGAL::AABB_traits Traits; + typedef CGAL::AABB_tree Tree; + typedef boost::optional::Type> Ray_intersection; + + cgal_shape_t polyhedron, polyhedron2, flattened; + + remove_thickness(const cgal_shape_t& p) + // edge_collapse(p) still does not work :( + : polyhedron(p) + , polyhedron2(p) + { + CGAL::Polygon_mesh_processing::triangulate_faces(polyhedron); + CGAL::Polygon_mesh_processing::triangulate_faces(polyhedron2); + + std::list non_degenerate, longitudonal; + std::set thin_sides; + + std::wcout << "ALL FACES:" << std::endl; + + for (auto& f : faces(polyhedron)) { + dump_facet(f); + if (facet_area(f) > 1.e-20) { + non_degenerate.push_back(f); + } + } + + std::wcout << "NON DEGENERATE:" << std::endl; + for (auto& f : non_degenerate) { + dump_facet(f); + } + + cgal_shape_t enlarged_non_degenerate_triangles; + Build_Offset bo; + bo.input = non_degenerate; + enlarged_non_degenerate_triangles.delegate(bo); + + Tree tree(non_degenerate.begin(), non_degenerate.end(), polyhedron); + + std::map face_normals; + boost::associative_property_map> face_normals_map(face_normals); + CGAL::Polygon_mesh_processing::compute_face_normals(polyhedron, face_normals_map); + + for (auto& f : non_degenerate) { + auto O = CGAL::centroid( + f->facet_begin()->vertex()->point(), + f->facet_begin()->next()->vertex()->point(), + f->facet_begin()->next()->next()->vertex()->point() + ); + + Ray ray(O, -face_normals_map[f]); + Ray_intersection intersection = tree.first_intersection(ray, [f](const cgal_shape_t::Facet_handle& p) { + return p == f; + }); + if (intersection) { + if (boost::get(&(intersection->first))) { + const Point* p = boost::get(&(intersection->first)); + const double d = std::sqrt(CGAL::to_double((*p - O).squared_length())); + if (d > 1.e-4) { + thin_sides.insert(f); + } + } + } else { + std::wcout << "No intersection :((!!!" << std::endl; + } + } + + std::wcout << "THIN SIDES:" << std::endl; + for (auto& f : thin_sides) { + dump_facet(f); + } + + for (auto& f : non_degenerate) { + if (thin_sides.find(f) == thin_sides.end()) { + longitudonal.push_back(f); + } + } + + std::wcout << "LONGITUDONAL:" << std::endl; + for (auto& f : longitudonal) { + dump_facet(f); + } + + cgal_shape_t enlarged_indiv_triangles; + Build_Offset bo; + bo.input = longitudonal; + enlarged_indiv_triangles.delegate(bo); + + { + std::ofstream ofs("enlarged.off"); + ofs.precision(17); + ofs << enlarged_indiv_triangles; + } + + Tree tree2(faces(enlarged_indiv_triangles).begin(), faces(enlarged_indiv_triangles).end(), enlarged_indiv_triangles); + + // std::map face_normals_2; + // boost::associative_property_map> face_normals_map_2(face_normals_2); + // CGAL::Polygon_mesh_processing::compute_face_normals(polyhedron, face_normals_map_2); + + // below does not seem to work? Do manually? + // std::map vertex_normals; + // boost::associative_property_map> vertex_normals_map(vertex_normals); + // CGAL::Polygon_mesh_processing::compute_normals(polyhedron, vertex_normals_map, face_normals_map_2); + + std::map new_points; + + for (Polyhedron::Facet_iterator fit = polyhedron.facets_begin(); + fit != polyhedron.facets_end(); + ++fit) { + if (CGAL::collinear( + fit->halfedge()->vertex()->point(), + fit->halfedge()->next()->vertex()->point(), + fit->halfedge()->opposite()->vertex()->point())) { + std::wcout << "degenerate triangle" << std::endl; + } + } + + /* + std::list vertices; + for (auto& f : non_degenerate) { + CGAL::Face_around_target_circulator it(f->halfedge(), polyhedron), end(it); + do { + vertices.push_back((*it)->halfedge()->vertex()); + ++it; + } while (it != end); + }*/ + + + for (auto& v : vertices(polyhedron)) { + auto O = v->point(); + + Kernel_::Vector_3 norm; + Kernel_::Vector_3 accum; + int count = 0; + CGAL::Face_around_target_circulator it(v->halfedge(), polyhedron), end(it); + do { + cgal_shape_t::Facet_handle fh = (*it)->halfedge()->facet(); + + auto jt = std::find(non_degenerate.begin(), non_degenerate.end(), fh); + std::wcout << "non degen: " << (jt != non_degenerate.end()) << std::endl; + auto kt = std::find(thin_sides.begin(), thin_sides.end(), fh); + std::wcout << "thin side: " << (kt != thin_sides.end()) << std::endl; + if (jt != non_degenerate.end() && kt == thin_sides.end()) { + // else degenerate, prevent div by zero, do not incorporate in vnorm. + // or else part of thin side + + auto p0 = (*it)->facet_begin()->vertex()->point(); + auto p1 = (*it)->facet_begin()->next()->vertex()->point(); + auto p2 = (*it)->facet_begin()->next()->next()->vertex()->point(); + + { + std::ostringstream oss; + oss.precision(8); + oss << "p0 " << p0.cartesian(0) << " " << p0.cartesian(1) << " " << p0.cartesian(2) << "\n"; + oss << "p1 " << p1.cartesian(0) << " " << p1.cartesian(1) << " " << p1.cartesian(2) << "\n"; + oss << "p2 " << p2.cartesian(0) << " " << p2.cartesian(1) << " " << p2.cartesian(2) << "\n"; + auto osss = oss.str(); + std::wcout << osss.c_str() << std::endl; + } + + auto fnorm = CGAL::cross_product(p0 - p1, p2 - p1); + fnorm /= std::sqrt(CGAL::to_double(fnorm.squared_length())); + + // const auto& fnorm = face_normals_map_2[*it]; + std::ostringstream oss; + oss.precision(8); + oss << fnorm.cartesian(0) << " " << fnorm.cartesian(1) << " " << fnorm.cartesian(2); + auto osss = oss.str(); + std::wcout << osss.c_str() << std::endl; + accum += fnorm; + + ++count; + } + + ++it; + } while (it != end); + + norm = accum / count; + std::wcout << "count " << count << std::endl; + + if (count == 0) { + continue; + } + + // v->vertex_begin(); + Ray ray(O, norm); + std::ostringstream oss; + oss.precision(8); + oss << O << " -> " << norm; + auto osss = oss.str(); + std::wcout << osss.c_str() << std::endl; + //// skip does not work anymore because we have offset the facets + // auto skip = [this, &v](const cgal_shape_t::Facet_handle& p) { + // CGAL::Face_around_target_circulator it(v->halfedge(), polyhedron), end(it); + // do { + // if ((*it)->facet_begin()->facet() == p) { + // return true; + // } + // } while (++it != end); + // return false; + // }; + std::list intersections; + tree2.all_intersections(ray, std::back_inserter(intersections)); + double N = std::numeric_limits::infinity(); + Point P; + if (intersections.size()) { + for (auto& intersection : intersections) { + if (boost::get(&(intersection->first))) { + const Point* p = boost::get(&(intersection->first)); + const double d = std::sqrt(CGAL::to_double((*p - O).squared_length())); + if (d < N) { + N = d; + P = *p; + } + std::wcout << "intersection @ " << d << std::endl; + } + } + std::wcout << "-----------" << std::endl; + + new_points[O] = P; + } else { + std::wcout << "no intersection :(" << std::endl; + } + } + + /* + for (auto& fi : thin_sides) { + auto f_circ = fi->facet_begin(); + polyhedron2.erase_facet(f_circ); + } + */ + + auto connected = connected_faces(*longitudonal.begin(), thin_sides); + + Builder_With_Map b2; + b2.input = connected; + b2.mapping = new_points; + + flattened.delegate(b2); + } + }; +} + void fix_spaceboundaries(IfcParse::IfcFile& f, bool no_progress, bool quiet, bool stderr_progress) { - typedef std::vector> > nefs_t; - typedef CGAL::Box_intersection_d::Box_with_handle_d Box; + typedef std::list> > nefs_t; + typedef CGAL::Box_intersection_d::Box_with_handle_d Box; + // typedef CGAL::Box_intersection_d::Box_d Box; + // std::map id_map; ifcopenshell::geometry::settings settings; settings.set(ifcopenshell::geometry::settings::USE_WORLD_COORDS, false); @@ -1295,38 +1870,31 @@ void fix_spaceboundaries(IfcParse::IfcFile& f, bool no_progress, bool quiet, boo // Apply transformation for (auto &vertex : vertices(s)) { vertex->point() = vertex->point().transform(trsf).transform(trsf2); - /* std::ostringstream ss; ss << vertex->point().cartesian(0); auto sss = ss.str(); std::wcout << sss.c_str() << std::endl; - */ } } - std::wcout << 1 << std::endl; CGAL::Nef_polyhedron_3 nef; auto c = convert_to_nef(s, nef); if (c != 0) { std::wcout << "Error " << c << std::endl; continue; } - std::wcout << 2 << std::endl; nef = CGAL::minkowski_sum_3(nef, cube); - std::wcout << 3 << std::endl; + std::wcout << "product: " << geom_object->product() << std::endl; nefs.push_back({ geom_object->product(), nef }); - std::wcout << 4 << std::endl; - - Kernel_::RT inf(std::numeric_limits::infinity()); - Kernel_::RT min[3] = { +inf, +inf, +inf }; - Kernel_::RT max[3] = { -inf, -inf, -inf }; - Box b(min, max, nefs.end() - 1); + Box b(&*(nefs.rbegin())); + // id_map[b.id()] = ; + for (auto &vertex : vertices(s)) { - Kernel_::RT p[3] = { - vertex->point().cartesian(0), - vertex->point().cartesian(1), - vertex->point().cartesian(2) + double p[3] = { + CGAL::to_double(vertex->point().cartesian(0)), + CGAL::to_double(vertex->point().cartesian(1)), + CGAL::to_double(vertex->point().cartesian(2)) }; b.extend(p); } @@ -1362,10 +1930,229 @@ void fix_spaceboundaries(IfcParse::IfcFile& f, bool no_progress, bool quiet, boo CGAL::box_self_intersection_d(boxes.begin(), boxes.end(), [](const Box& a, const Box& b) { std::ostringstream ss; - ss << a.handle()->first->data().toString() << "x" << b.handle()->first->data().toString() << std::endl; + // ss << id_map[a.id()]->first->data().toString() << "x" << id_map[b.id()]->first->data().toString() << std::endl; + // auto x = id_map[a.id()]->second * id_map[b.id()]->second; + + ss << a.handle()->first->data().toString() << "x" << a.handle()->first->data().toString() << std::endl; auto x = a.handle()->second * b.handle()->second; cgal_shape_t x_poly; x.convert_to_polyhedron(x_poly); + + CGAL::Polygon_mesh_processing::triangulate_faces(x_poly); + + auto vs = vertices(x_poly); + if (std::distance(vs.begin(), vs.end()) == 0) { + return; + } + + auto s0 = a.handle()->first->declaration().name(); + auto s1 = b.handle()->first->declaration().name(); + auto i0 = a.handle()->first->data().id(); + auto i1 = b.handle()->first->data().id(); + + if (s0 < s1) { + std::swap(s1, s0); + std::swap(i0, i1); + } + + { + auto FN = s0 + "-" + s1 + "-" + std::to_string(i0) + "-" + std::to_string(i1) + "sb.off"; + + std::ofstream os(FN.c_str()); + os.precision(17); + os << x_poly; + } + + remove_thickness r(x_poly); + + { + auto FN = s0 + "-" + s1 + "-" + std::to_string(i0) + "-" + std::to_string(i1) + "-sides-sb.off"; + + std::ofstream os(FN.c_str()); + os.precision(17); + os << r.polyhedron2; + } + + { + auto FN = s0 + "-" + s1 + "-" + std::to_string(i0) + "-" + std::to_string(i1) + "-flat-sb.off"; + + std::ofstream os(FN.c_str()); + os.precision(17); + os << r.flattened; + } + // r(); + + /* + std::map collapsed; + std::map collapsed_v; + + for (auto it = x_poly.edges_begin(); it != x_poly.edges_end(); ++it) { + auto& e = *it; + cgal_shape_t::Vertex_iterator v0 = e.vertex(); + cgal_shape_t::Vertex_iterator v1 = e.prev()->vertex(); + auto p0 = v0->point(); + auto p1 = v1->point(); + auto l = std::sqrt(CGAL::to_double((p1 - p0).squared_length())); + std::wcout << "edge w/ length " << l << std::endl; + + cgal_shape_t::Plane_3 plane(it->vertex()->point(), + it->next()->vertex()->point(), + it->next()->next()->vertex()->point()); + + auto d0 = plane.to_2d(e.prev()->vertex()->point()) - plane.to_2d(e.prev()->prev()->vertex()->point()); + auto d1 = plane.to_2d(e.vertex()->point()) - plane.to_2d(e.prev()->vertex()->point()); + auto d2 = plane.to_2d(e.next()->vertex()->point()) - plane.to_2d(e.vertex()->point()); + auto a0 = std::atan2(CGAL::to_double(d0.cartesian(1)), CGAL::to_double(d0.cartesian(0))); + auto a1 = std::atan2(CGAL::to_double(d1.cartesian(1)), CGAL::to_double(d1.cartesian(0))); + auto a2 = std::atan2(CGAL::to_double(d2.cartesian(1)), CGAL::to_double(d2.cartesian(0))); + auto a10 = a1 - a0; + auto a21 = a2 - a1; + if (a10 < 0.) { + a10 += 2 * M_PI; + } + if (a21 < 0.) { + a21 += 2 * M_PI; + } + + const bool is_convex = a10 < M_PI && a21 < M_PI; + + cgal_shape_t::Plane_3 opposite_plane(it->opposite()->vertex()->point(), + it->opposite()->next()->vertex()->point(), + it->opposite()->next()->next()->vertex()->point() + ); + + { + std::ostringstream oss; + oss << plane << " vs " << opposite_plane << "\n"; + oss << plane.orthogonal_vector() << " vs " << opposite_plane.orthogonal_vector(); + auto osss = oss.str(); + std::wcout << osss.c_str() << std::endl; + } + + bool is_internal = false; + if (std::sqrt(CGAL::to_double(plane.orthogonal_vector().squared_length())) < 1.e-15 || + std::sqrt(CGAL::to_double(opposite_plane.orthogonal_vector().squared_length())) < 1.e-15 + ) { + is_internal = true; + std::wcout << "Degenerate" << std::endl; + } else { + const double face_normal_dot = CGAL::to_double(approx_normalized(plane.orthogonal_vector()) * approx_normalized(opposite_plane.orthogonal_vector())); + std::wcout << "Face normal dot " << face_normal_dot << std::endl; + is_internal = face_normal_dot > 0.9; + } + + std::wcout << "Angles " << a0 << " " << a1 << " " << a2 << std::endl; + + if (l < 4.e-5 && (is_convex || is_internal)) { + auto p2 = CGAL::ORIGIN + ((p0 - CGAL::ORIGIN) + (p1 - CGAL::ORIGIN)) / 2; + + std::wcout << "(a) " << CGAL::to_double(p0.cartesian(0)) << " " << CGAL::to_double(p0.cartesian(1)) << " " << CGAL::to_double(p0.cartesian(2)) << "\n"; + std::wcout << "(b) " << CGAL::to_double(p1.cartesian(0)) << " " << CGAL::to_double(p1.cartesian(1)) << " " << CGAL::to_double(p1.cartesian(2)) << "\n"; + std::wcout << "(c) " << CGAL::to_double(p2.cartesian(0)) << " " << CGAL::to_double(p2.cartesian(1)) << " " << CGAL::to_double(p2.cartesian(2)) << "\n"; + // collapsed.insert({ v0, p2 }); + // collapsed.insert({ v1, p2 }); + collapsed.insert({ it, p2 }); + // Edges includes only half of the halfedges + collapsed.insert({ it->opposite(), p2 }); + + collapsed_v.insert({ v0, p2 }); + collapsed_v.insert({ v1, p2 }); + } + } + + { + auto FN = s0 + "-" + s1 + "-" + std::to_string(i0) + "-" + std::to_string(i1) + "sb.obj"; + std::ofstream ofs(FN.c_str()); + ofs.precision(17); + + int N = 1; + + std::set > faces_emitted; + for (auto& f : faces(x_poly)) { + std::ostringstream oss; + auto start = f->facet_begin(); + + bool part_collapsed = false; + + CGAL::Polyhedron_3::Halfedge_around_facet_const_circulator e = f->facet_begin(); + do { + auto it = collapsed.find(e); + if (it != collapsed.end()) { + part_collapsed = true; + break; + } + ++e; + } while (e != f->facet_begin()); + + decltype(faces_emitted)::key_type vss; + std::list points; + + if (!part_collapsed) { + e = f->facet_begin(); + do { + cgal_shape_t::Vertex_const_handle v = e->vertex(); + auto it = collapsed_v.find(v); + if (it == collapsed_v.end()) { + std::wcout << "Unexpected " + << CGAL::to_double(v->point().cartesian(0)) << " " + << CGAL::to_double(v->point().cartesian(1)) << " " + << CGAL::to_double(v->point().cartesian(2)) << std::endl; + } else { + points.push_back(it->second); + vss.insert(it->second); + } + ++e; + } while (e != f->facet_begin()); + + if (faces_emitted.find(vss) != faces_emitted.end()) { + std::wcout << "Emitted" << std::endl; + } else { + faces_emitted.insert(vss); + + for (auto& p : points) { + ofs << "v " << CGAL::to_double(p.cartesian(0)) << " " << CGAL::to_double(p.cartesian(1)) << " " << CGAL::to_double(p.cartesian(2)) << "\n"; + } + + ofs << "f "; + for (auto i = 0; i < points.size(); ++i) { + if (i) { + ofs << " "; + } + ofs << i + N; + } + ofs << "\n"; + + N += points.size(); + } + } + } + } + + */ + + /* + // edge collapse does not work on the rational number types + typedef CGAL::Simple_cartesian simple; + CGAL::Polyhedron_3 x_simple; + poly_copy(x_simple, x_poly); + + // flattening from a thin box to a plane is not valid in edge_collapse() + Stats stats; + My_visitor vis(&stats); + SMS::Edge_length_cost elc; + SMS::Edge_length_stop_predicate stop(1.e-3); + + int r = SMS::edge_collapse(x_simple, stop, + CGAL::parameters::vertex_index_map(get(CGAL::vertex_external_index, x_simple)) + .halfedge_index_map(get(CGAL::halfedge_external_index, x_simple)) + .visitor(vis) + .get_cost(elc) + ); + + std::wcout << "Removed: " << r << std::endl; + */ + + /* for (auto& v : vertices(x_poly)) { auto p = v->point(); for (int i = 0; i < 3; ++i) { @@ -1376,6 +2163,7 @@ void fix_spaceboundaries(IfcParse::IfcFile& f, bool no_progress, bool quiet, boo ss << "---" << std::endl; auto sss = ss.str(); std::wcout << sss.c_str(); + */ }); if (!no_progress && quiet) { diff --git a/src/ifcgeom/kernel_agnostic/blog.txt b/src/ifcgeom/kernel_agnostic/blog.txt deleted file mode 100644 index e188c8d407..0000000000 --- a/src/ifcgeom/kernel_agnostic/blog.txt +++ /dev/null @@ -1,71 +0,0 @@ -v0.6.0 - -People not following the development of IfcOpenShell actively and happily using the master branch of the github repository might be surprised to know there is a lot of activity happening in the v0.6.0 and v0.7.0 branches. This post discusses the changes in the v0.6.0 branch. The following post will elaborate on some of the design decisions we are making in the v0.7.0 branch. - -Schemas - -The most significant improvement in the v0.6.0 branch is that multiple schemas (IFC2X3, IFC4, IFC4X1 and IFC4X2) are supported from within the same executable, module or plug-in. Previously, selecting the schema had been a compile-time option. - -In IfcOpenShell and most other EXPRESS-based toolkits, the IFC schema is compiled into (a) the early-bound definitions: a class hierarchy with member functions and (b) a set of methods to operate on the schema definitions at runtime (late-bound access). C++ only allows very limited introspection (but the development of C++ is very active, see for example P1240 https://github.com/cplusplus/papers/issues/545) so to complement the lack of introspection a set of methods exists to query for example all attribute names or the sub- and supertypes of an entity. In the master branch these methods are static, in the v0.6.0 branch these are the member functions of a schema class, that is a more complete reference mirrorring the EXPRESS schema definition at runtime. See IfcBaseEntity::declararation() or IfcParse::schema::declaration_by_name("IfcWall")->as_entity()->all_attribute_names(). - -Writing schema agnostic code - -The code generated from the four schemas are completely orthogonal class hiercharies. For the C++ compiler there is no relationship between a Ifc2x3::IfcWall and a Ifc4::IfcWall. But IfcOpenShell offers three ways to write code that adapts to the schema of the file known at runtime. - -(a) preprocessor - -This is the approach taken in the IfcGeom modules in v0.6.0. Essentially the same code base is compiled multiple times where the schema is available as a preprocessor constant. This means you can enable specific code paths with for example #ifdef directives. In this way the added entities in Ifc4 (IfcBSplineSurface, yay!) can be selectively compiled for example. - -https://github.com/IfcOpenShell/IfcOpenShell/blob/v0.6.0/src/ifcgeom/IfcGeomFaces.cpp#L1127 - -Smaller code blocks can be written as macros as well. - -https://github.com/IfcOpenShell/IfcOpenShell/blob/v0.6.0/src/ifcgeom_schema_agnostic/Kernel.cpp#L74 - -Benefits: fairly readible code, full autocompletion typically in an IDE when using the static library approach -Downsides: Some infrastructure required to compile the different libraries and select the correct implementation at runtime - -(b) late-bound access - -There are two modes of accessing schemas. In the early-bound approach function signatures and return types are known at compilation time. In the late-bound approach attribute names are referenced by strings and types are - -Ifc2x3::IfcWall* wall; -// Early-bound access; -std::string global_id = wall->GlobalId(); -// Late-bound access. -std::string global_id = *wall->get("GlobalId"); -// ERROR: By dereferencing the return type, it is casted into a string, which will cause an exception *at runtime* when the types do not match. -int global_id = *wall->get("GlobalId"); - -Benefits: -fairly readible code -no complicated setup of different libraries -Downsides: -no code completion -errors are only spotted at runtime, not compile-time -late-bound manipulation of inverse attributes is not well supported currently in IfcOpenShell -less means for the compiler to create highly optimized code - -(c) templates - -C++ has very extensive support for compile time generic arguments: templates. - -template -void print_globalid(Schema::IfcWall* wall) { - std::cout << wall->GlobalId(); -} - -Benefits: -no complicated setup of different libraries -no autocompletion typically, but errors caught at compile-time -Downsides: -fairly unreadible code due to additional template and typename keywords. -error messages are harder to make sense up (due to two phase lookup rules for example) - -All three approaches are used in the IfcOpenShell code-base. - -Other improvements: - -Multi-threading in collaboration with TNO, MAUC and Airsquire - -Direct binary glTF output (previously supported through Collada and Collada2Gltf) in collaboration with Schuco US. diff --git a/src/ifcgeom/schema/mapping.cpp b/src/ifcgeom/schema/mapping.cpp index cbc395c51d..1e063ec52d 100644 --- a/src/ifcgeom/schema/mapping.cpp +++ b/src/ifcgeom/schema/mapping.cpp @@ -700,12 +700,15 @@ void mapping::get_representations(std::vector& tasks, continue; } - geometry_conversion_task task; - task.index = task_index++; - task.representation = representation; - task.products = ifcproducts->generalize(); + // @todo, fix this properly by considering the mapped geometry types in the representation. + if (representation->hasRepresentationIdentifier() && representation->RepresentationIdentifier() == "Body") { + geometry_conversion_task task; + task.index = task_index++; + task.representation = representation; + task.products = ifcproducts->generalize(); - tasks.emplace_back(task); + tasks.emplace_back(task); + } } }