From 4374c81981567c37223731e5baafaf69794bca07 Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Tue, 1 Sep 2026 10:35:11 -0700 Subject: [PATCH] Honour IfcAxis2PlacementLinear Axis/RefDirection in the loft builder (#9408) * Honour IfcAxis2PlacementLinear Axis/RefDirection in the loft builder make_loft() (src/ifcgeom/infra_sweep_helper.cpp), shared by IfcSectionedSolidHorizontal and IfcSectionedSurface, mishandled a cross section's IfcAxis2PlacementLinear in two ways: 1. A placement carrying Axis but no RefDirection was placed with a fixed [e_y | e_z | e_x] world-axis permutation that ignored the directrix. On any directrix not running along +X (e.g. a north-south road pavement, or anywhere along a curve) the profile came out mis-oriented or collapsed to a sliver. 2. When two adjacent CrossSectionPositions used direction vectors inconsistently (a raked RefDirection at one, a plain Axis at the other) make_loft() logged GEO 42, dropped the rotation for the whole segment and squared every cap -- and in one configuration left the sweep frame flipped, so OpenCASCADE failed to build the solid at all. Now a small profile_basis() helper builds every cross section's frame the same way: profile Y = Axis, profile normal = RefDirection, and -- when RefDirection is absent -- profile normal = the directrix tangent, so the section stays perpendicular to the path (buildingSMART IFC4.x-IF #147). When the two bracketing placements ask for the same orientation the sweep frame carries it, built against the curve. When they disagree the sweep frame stays on the shared Axis (continuous with the neighbouring consistent segments, so nothing flips) and each end's own authored orientation is folded into its profile points via a change of basis, so each end cap still lands exactly as authored while the body in between keeps following the directrix. The all-equal and no-direction-vector paths are unchanged. The two mappings now also carry the raw RefDirection through on cross_section, alongside the existing rotation matrix. Adds C++ tests (a raked end logs no GEO 42; a directrix that does not run along +X still lofts a full-size solid) and Python tests (uniform prism raked at one end and square at the other; a north-south directrix keeps its width; OffsetLateral/OffsetVertical are scaled by the model length unit). * Renames profile_rotations to profile_axis for consistency with profile_ref_directions --- src/ifcgeom/infra_sweep_helper.cpp | 120 ++++++-- src/ifcgeom/infra_sweep_helper.h | 6 + .../mapping/IfcSectionedSolidHorizontal.cpp | 25 +- src/ifcgeom/mapping/IfcSectionedSurface.cpp | 31 +- .../tests/test_ifcopenshell_geometry.cpp | 146 ++++++++++ .../test/test_create_shape.py | 272 ++++++++++++++++++ 6 files changed, 558 insertions(+), 42 deletions(-) diff --git a/src/ifcgeom/infra_sweep_helper.cpp b/src/ifcgeom/infra_sweep_helper.cpp index 3c65dd3d60..b6c0dc109b 100644 --- a/src/ifcgeom/infra_sweep_helper.cpp +++ b/src/ifcgeom/infra_sweep_helper.cpp @@ -2,7 +2,9 @@ #include "infra_sweep_helper.h" #include "function_item_evaluator.h" +#include #include +#include using namespace ifcopenshell::geom; @@ -72,6 +74,53 @@ taxonomy::loft::ptr ifcopenshell::geom::make_loft(const ifcopenshell::geom::sett longitudes.push_back(x.dist_along); } longitudes.push_back(std::numeric_limits::infinity()); + + // Directrix frame (col0 = tangent, col1 = lateral, col2 = up, col3 = position) at + // every cross section station. Only needed to reconcile a cross section's own + // Axis / RefDirection against the curve; skip the work when no placement has one. + std::vector section_directrix_frames; + if (std::any_of(cross_sections.begin(), cross_sections.end(), + [](const cross_section& cs) { return cs.rotation.has_value(); })) { + section_directrix_frames.reserve(cross_sections.size()); + for (const auto& cs : cross_sections) { + section_directrix_frames.push_back(evaluator.evaluate(std::min(std::max(cs.dist_along, start), end))); + } + } + + // The frame a cross section is placed in, given the directrix frame at its station + // and its own IfcAxis2PlacementLinear: profile X, profile Y = Axis, profile normal + // = RefDirection -- or, when RefDirection was not authored, the curve tangent, so + // the section stays perpendicular to the path (buildingSMART IFC4.x-IF #147). When + // the placement carries no direction vectors the section just follows the curve + // (lateral, up, tangent). + const auto profile_basis = + [](const std::optional& rotation, + const std::optional& ref_direction, + const Eigen::Matrix4d& directrix_frame) -> Eigen::Matrix3d { + const Eigen::Vector3d tangent = directrix_frame.col(0).head<3>().normalized(); + const Eigen::Vector3d lateral = directrix_frame.col(1).head<3>().normalized(); + const Eigen::Vector3d up = directrix_frame.col(2).head<3>().normalized(); + Eigen::Matrix3d B; + if (!rotation) { + B.col(0) = lateral; + B.col(1) = up; + B.col(2) = tangent; + return B; + } + const Eigen::Vector3d axis = rotation->col(2).normalized(); + const Eigen::Vector3d normal = ref_direction ? ref_direction->normalized() : tangent; + Eigen::Vector3d x = axis.cross(normal); + if (x.norm() < 1.e-9) { + // Axis parallel to the normal: fall back to the curve's own lateral. + x = lateral - lateral.dot(axis) * axis; + } + x.normalize(); + B.col(0) = x; + B.col(1) = axis; + B.col(2) = x.cross(axis); + return B; + }; + auto profile_index = longitudes.begin(); for (size_t i = 0; i <= num_steps; ++i) { auto dist_along = start + delta_step * i; @@ -88,6 +137,7 @@ taxonomy::loft::ptr ifcopenshell::geom::make_loft(const ifcopenshell::geom::sett const auto& profile_a = cross_sections[std::distance(longitudes.begin(), profile_index)].section_geometry; const auto& offset_a = cross_sections[std::distance(longitudes.begin(), profile_index)].offset; const auto& rotation_a = cross_sections[std::distance(longitudes.begin(), profile_index)].rotation; + const auto& ref_direction_a = cross_sections[std::distance(longitudes.begin(), profile_index)].ref_direction; taxonomy::geom_item::ptr interpolated = nullptr; @@ -98,26 +148,37 @@ taxonomy::loft::ptr ifcopenshell::geom::make_loft(const ifcopenshell::geom::sett (profile_index + 1 < longitudes.end()) && (relative_dist_along >= 1.e-9 || offset_a.cwiseAbs().maxCoeff() > 0. || rotation_a); + // When both bracketing placements ask for the same orientation, drive the sweep + // frame from it directly (relative to the curve). When they disagree, keep the + // sweep frame on the curve and fold each section's own orientation into its + // profile points via section_basis_a / section_basis_b so the ends still land + // exactly as authored without twisting the body between them. std::optional interpolated_rotation; + std::optional interpolated_ref_direction; + Eigen::Matrix3d section_basis_a = Eigen::Matrix3d::Identity(); + Eigen::Matrix3d section_basis_b = Eigen::Matrix3d::Identity(); if (should_interpolate) { taxonomy::geom_item::ptr profile_b; Eigen::Vector3d offset_b; std::optional rotation_b; + std::optional ref_direction_b; if ((profile_index + 1 < longitudes.end())) { profile_b = cross_sections[std::distance(longitudes.begin(), profile_index) + 1].section_geometry; offset_b = cross_sections[std::distance(longitudes.begin(), profile_index) + 1].offset; rotation_b = cross_sections[std::distance(longitudes.begin(), profile_index) + 1].rotation; + ref_direction_b = cross_sections[std::distance(longitudes.begin(), profile_index) + 1].ref_direction; } else { profile_b = profile_a; offset_b = offset_a; rotation_b = rotation_a; + ref_direction_b = ref_direction_a; } // Only interpolate if the profiles are different or either of the offsets is non-zero bool should_interpolate2 = (profile_a->instance != profile_b->instance) || - (offset_a.cwiseAbs().maxCoeff() > 0. || offset_b.cwiseAbs().maxCoeff() > 0. || rotation_b); + (offset_a.cwiseAbs().maxCoeff() > 0. || offset_b.cwiseAbs().maxCoeff() > 0. || rotation_a || rotation_b); if (should_interpolate2) { @@ -160,12 +221,33 @@ taxonomy::loft::ptr ifcopenshell::geom::make_loft(const ifcopenshell::geom::sett } auto interpolated_offset = lerp(offset_a, offset_b, relative_dist_along); - if (rotation_a == rotation_b && rotation_a) { - // @todo we don't support an overridden rotation on only one of the placements - // in which case we would need to lerp with the rotation component below in m4b. - interpolated_rotation = lerp(*rotation_a, *rotation_b, relative_dist_along); - } else if (rotation_a != rotation_b) { - logger.error("GEO", 42, "Direction vectors on cross section placements only supported when used consistently"); + if (rotation_a == rotation_b) { + // Same orientation on both placements (including both absent): the + // sweep frame carries it, built against the curve just below. + interpolated_rotation = rotation_a; + interpolated_ref_direction = ref_direction_a; + } else if (rotation_a || rotation_b) { + // The two placements disagree -- in practice they share an Axis but + // only one carries a RefDirection (a raked end against a square run). + // Drive the sweep frame from the shared Axis with the profile normal + // on the curve tangent -- identical to the neighbouring consistent + // segments, so m4b stays continuous across the boundary and the body + // never flips -- then fold each section's *own* authored orientation + // into its profile points through section_basis_a / section_basis_b, + // so each end cap still lands exactly as authored. + const auto ia = static_cast(std::distance(longitudes.begin(), profile_index)); + assert(ia + 1 < section_directrix_frames.size()); + const std::optional no_ref; + const auto base_a = profile_basis(rotation_a, no_ref, section_directrix_frames[ia]); + const auto base_b = profile_basis(rotation_b, no_ref, section_directrix_frames[ia + 1]); + section_basis_a = + base_a.transpose() * + profile_basis(rotation_a, ref_direction_a, section_directrix_frames[ia]); + section_basis_b = + base_b.transpose() * + profile_basis(rotation_b, ref_direction_b, section_directrix_frames[ia + 1]); + interpolated_rotation = rotation_a ? rotation_a : rotation_b; + interpolated_ref_direction = no_ref; } taxonomy::loop::ptr w1, w2; @@ -337,7 +419,9 @@ taxonomy::loft::ptr ifcopenshell::geom::make_loft(const ifcopenshell::geom::sett const auto& tagged_point_on_w1 = tag_to_point_on_w1[t]; const auto& tagged_point_on_w2 = tag_to_point_on_w2[t]; - auto p3 = (lerp(tagged_point_on_w1->ccomponents(), tagged_point_on_w2->ccomponents(), relative_dist_along) + interpolated_offset).eval(); + const Eigen::Vector3d rebased_w1 = section_basis_a * tagged_point_on_w1->ccomponents(); + const Eigen::Vector3d rebased_w2 = section_basis_b * tagged_point_on_w2->ccomponents(); + auto p3 = (lerp(rebased_w1, rebased_w2, relative_dist_along) + interpolated_offset).eval(); std::set tags_for_this_point_on_subsequent_profile = {t}; @@ -357,7 +441,9 @@ taxonomy::loft::ptr ifcopenshell::geom::make_loft(const ifcopenshell::geom::sett } else { for (auto tmp__ : boost::combine(w1_points, w2_points)) { boost::tie(p1, p2) = tmp__; - auto p3 = (lerp(p1->ccomponents(), p2->ccomponents(), relative_dist_along) + interpolated_offset).eval(); + const Eigen::Vector3d rebased_1 = section_basis_a * p1->ccomponents(); + const Eigen::Vector3d rebased_2 = section_basis_b * p2->ccomponents(); + auto p3 = (lerp(rebased_1, rebased_2, relative_dist_along) + interpolated_offset).eval(); points.push_back(taxonomy::make(p3)); } } @@ -399,17 +485,13 @@ taxonomy::loft::ptr ifcopenshell::geom::make_loft(const ifcopenshell::geom::sett std::wcout << "#" << pwf->instance.data().id() << " " << dist_along << ": " << m4.col(3).row(2).value() << std::endl; }*/ + // Sweep frame at this station: the profile orientation asked for by the + // (consistent) placements, built against the curve here so it follows the + // directrix. Falls back to the plain curve frame (lateral, up, tangent) when + // no placement carries direction vectors. Inconsistent placements keep this + // on the curve and are reconciled through section_basis_a / section_basis_b. Eigen::Matrix4d m4b = Eigen::Matrix4d::Identity(); - if (interpolated_rotation) { - // direction vectors on the linear placement overwrite the placement otherwise inferred from the tangent - m4b.col(0).head<3>() = interpolated_rotation->col(1); - m4b.col(1).head<3>() = interpolated_rotation->col(2); - m4b.col(2).head<3>() = interpolated_rotation->col(0); - } else { - m4b.col(0).head<3>() = m4.col(1).head<3>().normalized(); - m4b.col(1).head<3>() = m4.col(2).head<3>().normalized(); - m4b.col(2).head<3>() = m4.col(0).head<3>().normalized(); - } + m4b.block<3, 3>(0, 0) = profile_basis(interpolated_rotation, interpolated_ref_direction, m4); m4b.col(3).head<3>() = m4.col(3).head<3>(); if (interpolated) { diff --git a/src/ifcgeom/infra_sweep_helper.h b/src/ifcgeom/infra_sweep_helper.h index cb0c06ac39..709292a2ca 100644 --- a/src/ifcgeom/infra_sweep_helper.h +++ b/src/ifcgeom/infra_sweep_helper.h @@ -14,7 +14,13 @@ namespace ifcopenshell { double dist_along; taxonomy::geom_item::ptr section_geometry; Eigen::Vector3d offset; + // rotation: the IfcAxis2PlacementLinear basis [X | Y | Axis], or nullopt when + // the placement carries neither Axis nor RefDirection. std::optional rotation; + // ref_direction: the raw RefDirection, only when it was actually authored on the + // placement (rotation carries a default otherwise). When absent the profile + // normal comes from the directrix tangent -- see buildingSMART IFC4.x-IF #147. + std::optional ref_direction; bool operator <(const cross_section& other) const { return dist_along < other.dist_along; diff --git a/src/ifcgeom/mapping/IfcSectionedSolidHorizontal.cpp b/src/ifcgeom/mapping/IfcSectionedSolidHorizontal.cpp index 546c8712d6..e790b36a35 100644 --- a/src/ifcgeom/mapping/IfcSectionedSolidHorizontal.cpp +++ b/src/ifcgeom/mapping/IfcSectionedSolidHorizontal.cpp @@ -49,7 +49,8 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcSectionedSolidHorizontal& in // The longitudes determine the range of the sweep and the offsets are interpolated in between // sweep segments. std::vector profile_offsets; - std::vector> profile_rotations; + std::vector> profile_axes; + std::vector> profile_ref_directions; std::vector longitudes; for (auto& cs : css) { @@ -72,23 +73,27 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcSectionedSolidHorizontal& in profile_offsets.push_back(po); std::optional rot; + std::optional ref_direction; if (csp.Axis() && csp.RefDirection()) { + ref_direction = taxonomy::cast(map(csp.RefDirection()))->ccomponents(); rot = taxonomy::matrix4( Eigen::Vector3d(0, 0, 0), taxonomy::cast(map(csp.Axis()))->ccomponents(), - taxonomy::cast(map(csp.RefDirection()))->ccomponents()).ccomponents().block<3,3>(0,0); + *ref_direction).ccomponents().block<3,3>(0,0); } else if (csp.Axis()) { rot = taxonomy::matrix4( Eigen::Vector3d(0, 0, 0), taxonomy::cast(map(csp.Axis()))->ccomponents()).ccomponents().block<3, 3>(0, 0); - } else if (csp.RefDirection()) { - rot = taxonomy::matrix4( - Eigen::Vector3d(0, 0, 0), - Eigen::Vector3d(0, 0, 1), - taxonomy::cast(map(csp.RefDirection()))->ccomponents() + } else if (csp.RefDirection()) { + ref_direction = taxonomy::cast(map(csp.RefDirection()))->ccomponents(); + rot = taxonomy::matrix4( + Eigen::Vector3d(0, 0, 0), + Eigen::Vector3d(0, 0, 1), + *ref_direction ).ccomponents().block<3, 3>(0, 0); - } - profile_rotations.push_back(rot); + } + profile_axes.push_back(rot); + profile_ref_directions.push_back(ref_direction); } if (faces.size() != profile_offsets.size()) { logger_.warning("GEO", 286, "Expected CrossSections and CrossSectionPositions to be equal length, but got " + std::to_string(faces.size()) + " and " + std::to_string(profile_offsets.size()) + " respectively", inst); @@ -100,7 +105,7 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcSectionedSolidHorizontal& in } for (size_t i = 0; i < faces.size(); ++i) { - cross_sections.push_back({ longitudes[i], faces[i], profile_offsets[i], profile_rotations[i]}); + cross_sections.push_back({longitudes[i], faces[i], profile_offsets[i], profile_axes[i], profile_ref_directions[i]}); } #else return nullptr; diff --git a/src/ifcgeom/mapping/IfcSectionedSurface.cpp b/src/ifcgeom/mapping/IfcSectionedSurface.cpp index a0ef5f727b..32423d9f7e 100644 --- a/src/ifcgeom/mapping/IfcSectionedSurface.cpp +++ b/src/ifcgeom/mapping/IfcSectionedSurface.cpp @@ -50,7 +50,8 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcSectionedSurface& inst) { // The longitudes determine the range of the sweep and the offsets are interpolated in between // sweep segments. std::vector profile_offsets; - std::vector> profile_rotations; + std::vector> profile_axes; + std::vector> profile_ref_directions; std::vector longitudes; for (auto& cs : css) { @@ -74,24 +75,28 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcSectionedSurface& inst) { profile_offsets.push_back(po); std::optional rot; + std::optional ref_direction; if (csp.Axis() && csp.RefDirection()) { + ref_direction = taxonomy::cast(map(csp.RefDirection()))->ccomponents(); rot = taxonomy::matrix4( Eigen::Vector3d(0, 0, 0), taxonomy::cast(map(csp.Axis()))->ccomponents(), - taxonomy::cast(map(csp.RefDirection()))->ccomponents()).ccomponents().block<3, 3>(0, 0); + *ref_direction).ccomponents().block<3, 3>(0, 0); } else if (csp.Axis()) { rot = taxonomy::matrix4( Eigen::Vector3d(0, 0, 0), taxonomy::cast(map(csp.Axis()))->ccomponents()).ccomponents().block<3, 3>(0, 0); - } else if (csp.RefDirection()) { - rot = taxonomy::matrix4( - Eigen::Vector3d(0, 0, 0), - Eigen::Vector3d(0, 0, 1), - taxonomy::cast(map(csp.RefDirection()))->ccomponents()) - .ccomponents() - .block<3, 3>(0, 0); - } - profile_rotations.push_back(rot); + } else if (csp.RefDirection()) { + ref_direction = taxonomy::cast(map(csp.RefDirection()))->ccomponents(); + rot = taxonomy::matrix4( + Eigen::Vector3d(0, 0, 0), + Eigen::Vector3d(0, 0, 1), + *ref_direction) + .ccomponents() + .block<3, 3>(0, 0); + } + profile_axes.push_back(rot); + profile_ref_directions.push_back(ref_direction); } #else return nullptr; @@ -106,8 +111,8 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcSectionedSurface& inst) { } for (size_t i = 0; i < faces.size(); ++i) { - cross_sections.push_back({ longitudes[i], faces[i], profile_offsets[i], profile_rotations[i] }); - } + cross_sections.push_back({longitudes[i], faces[i], profile_offsets[i], profile_axes[i], profile_ref_directions[i]}); + } } return make_loft(settings_, inst, fn, cross_sections); diff --git a/src/ifcgeom/tests/test_ifcopenshell_geometry.cpp b/src/ifcgeom/tests/test_ifcopenshell_geometry.cpp index 63f20538dc..7933bcb8cf 100644 --- a/src/ifcgeom/tests/test_ifcopenshell_geometry.cpp +++ b/src/ifcgeom/tests/test_ifcopenshell_geometry.cpp @@ -1,7 +1,9 @@ #include +#include #include #include +#include #include #include "ifcgeom/converter.h" @@ -148,8 +150,152 @@ const ifcopenshell::geom::geometry_conversion_task* task_for_product( return nullptr; } +// A minimal IfcSectionedSolidHorizontal (IFC4X3_ADD2) whose two +// IfcAxis2PlacementLinear cross section positions use direction vectors +// inconsistently: the near position carries a raked RefDirection and a +// 1/cos(theta) wider profile, the far position carries neither. Before the +// make_loft() fix this logged GEO 42 and dropped the rotation, lofting a wedge. +constexpr const char* RAKED_SECTIONED_SOLID_SPF = R"IFC(ISO-10303-21; +HEADER; +FILE_DESCRIPTION((''),'2;1'); +FILE_NAME('','',(''),(''),'','',''); +FILE_SCHEMA(('IFC4X3_ADD2')); +ENDSEC; +DATA; +#1=IFCPROJECT('0RYK8PV8D0ee9DDm77xcTZ',$,'T',$,$,$,$,(#6),$); +#2=IFCCARTESIANPOINT((0.,0.,0.)); +#3=IFCDIRECTION((0.,0.,1.)); +#4=IFCDIRECTION((1.,0.,0.)); +#5=IFCAXIS2PLACEMENT3D(#2,#3,#4); +#6=IFCGEOMETRICREPRESENTATIONCONTEXT($,'Model',3,1.E-05,#5,$); +#7=IFCGEOMETRICREPRESENTATIONSUBCONTEXT('Body','Model',*,*,*,*,#6,$,.MODEL_VIEW.,$); +#8=IFCCARTESIANPOINT((0.,0.,0.)); +#9=IFCCARTESIANPOINT((40.,0.,0.)); +#10=IFCPOLYLINE((#8,#9)); +#11=IFCDIRECTION((0.,0.,1.)); +#12=IFCDIRECTION((0.9034641832977311,-0.42866358545853134,0.)); +#13=IFCPOINTBYDISTANCEEXPRESSION(IFCLENGTHMEASURE(0.),$,$,$,#10); +#14=IFCAXIS2PLACEMENTLINEAR(#13,#11,#12); +#15=IFCPOINTBYDISTANCEEXPRESSION(IFCLENGTHMEASURE(40.),$,$,$,#10); +#16=IFCAXIS2PLACEMENTLINEAR(#15,#11,$); +#17=IFCCARTESIANPOINTLIST2D(((-0.9223756168081689,0.),(0.9223756168081689,0.),(0.9223756168081689,6.),(-0.9223756168081689,6.),(-0.9223756168081689,0.)),$); +#18=IFCINDEXEDPOLYCURVE(#17,$,.F.); +#19=IFCARBITRARYCLOSEDPROFILEDEF(.AREA.,$,#18); +#20=IFCCARTESIANPOINTLIST2D(((-0.8333333333333334,0.),(0.8333333333333334,0.),(0.8333333333333334,6.),(-0.8333333333333334,6.),(-0.8333333333333334,0.)),$); +#21=IFCINDEXEDPOLYCURVE(#20,$,.F.); +#22=IFCARBITRARYCLOSEDPROFILEDEF(.AREA.,$,#21); +#23=IFCSECTIONEDSOLIDHORIZONTAL(#10,(#19,#22),(#14,#16)); +#24=IFCBUILDINGELEMENTPROXY('3cWcr4$892GAeKwryhDILR',$,'wingwall',$,$,$,#26,$,$); +#25=IFCSHAPEREPRESENTATION(#7,'Body','AdvancedSweptSolid',(#23)); +#26=IFCPRODUCTDEFINITIONSHAPE($,$,(#25)); +ENDSEC; +END-ISO-10303-21; +)IFC"; + +struct sectioned_solid_result { + std::size_t geo42_count = 0; + bool produced_brep = false; + double projected_area_x = 0.0; + double projected_area_y = 0.0; + double projected_area_z = 0.0; +}; + +sectioned_solid_result convert_sectioned_solid(const std::string& spf) { + std::istringstream stream(spf); + ifcopenshell::logger log; + log.output_format(ifcopenshell::logger::FMT_INMEMORY); + ifcopenshell::file file(stream, static_cast(spf.size()), log); + REQUIRE(file.good()); + + ifcopenshell::geom::settings settings; + ifcopenshell::geom::converter converter( + ifcopenshell::geom::kernels::construct(&file, "opencascade", settings, log), &file, settings, log); + + std::vector tasks; + std::vector filters; + converter.mapping()->get_representations(tasks, filters); + REQUIRE(!tasks.empty()); + + sectioned_solid_result result; + for (const auto& task : tasks) { + REQUIRE(!task.products.empty()); + auto* elem = converter.create_brep_for_representation_and_product(task.representation, task.products.front()); + if (elem) { + result.produced_brep = true; + elem->calculate_projected_surface_area( + result.projected_area_x, result.projected_area_y, result.projected_area_z); + } + delete elem; + } + result.geo42_count = log.count("GEO42"); + return result; +} + +std::size_t count_geo42_converting(const std::string& spf) { + return convert_sectioned_solid(spf).geo42_count; +} + +// A minimal IfcSectionedSolidHorizontal (IFC4X3_ADD2) whose cross section +// placements carry an explicit Axis = (0,0,1) but no RefDirection, on a +// directrix that runs along +Y (not the global +X). Per buildingSMART +// IFC4.x-IF #147 the profile normal follows the directrix tangent, so the +// 12 x 0.5 rectangle sweeps 60 along +Y: a plan (Z) projected area of ~720. +// Before the fix the profile was placed with a fixed axis permutation that +// ignored the directrix and the solid collapsed. +constexpr const char* AXIS_ALIGNED_SECTIONED_SOLID_SPF = R"IFC(ISO-10303-21; +HEADER; +FILE_DESCRIPTION((''),'2;1'); +FILE_NAME('','',(''),(''),'','',''); +FILE_SCHEMA(('IFC4X3_ADD2')); +ENDSEC; +DATA; +#1=IFCPROJECT('0RYK8PV8D0ee9DDm77xcTZ',$,'T',$,$,$,$,(#6),$); +#2=IFCCARTESIANPOINT((0.,0.,0.)); +#3=IFCDIRECTION((0.,0.,1.)); +#4=IFCDIRECTION((1.,0.,0.)); +#5=IFCAXIS2PLACEMENT3D(#2,#3,#4); +#6=IFCGEOMETRICREPRESENTATIONCONTEXT($,'Model',3,1.E-05,#5,$); +#7=IFCGEOMETRICREPRESENTATIONSUBCONTEXT('Body','Model',*,*,*,*,#6,$,.MODEL_VIEW.,$); +#8=IFCCARTESIANPOINT((0.,0.,0.)); +#9=IFCCARTESIANPOINT((0.,60.,0.)); +#10=IFCPOLYLINE((#8,#9)); +#11=IFCDIRECTION((0.,0.,1.)); +#13=IFCPOINTBYDISTANCEEXPRESSION(IFCLENGTHMEASURE(0.),$,$,$,#10); +#14=IFCAXIS2PLACEMENTLINEAR(#13,#11,$); +#15=IFCPOINTBYDISTANCEEXPRESSION(IFCLENGTHMEASURE(60.),$,$,$,#10); +#16=IFCAXIS2PLACEMENTLINEAR(#15,#11,$); +#17=IFCCARTESIANPOINTLIST2D(((-6.,0.),(6.,0.),(6.,0.5),(-6.,0.5),(-6.,0.)),$); +#18=IFCINDEXEDPOLYCURVE(#17,$,.F.); +#19=IFCARBITRARYCLOSEDPROFILEDEF(.AREA.,$,#18); +#23=IFCSECTIONEDSOLIDHORIZONTAL(#10,(#19,#19),(#14,#16)); +#24=IFCBUILDINGELEMENTPROXY('3cWcr4$892GAeKwryhDILR',$,'pavement',$,$,$,#26,$,$); +#25=IFCSHAPEREPRESENTATION(#7,'Body','AdvancedSweptSolid',(#23)); +#26=IFCPRODUCTDEFINITIONSHAPE($,$,(#25)); +ENDSEC; +END-ISO-10303-21; +)IFC"; + } // namespace +TEST_CASE("IfcSectionedSolidHorizontal raked end cut does not log GEO 42", "[ifcgeom][infra-sweep]") { + if (std::string(STRINGIFY(IfcSchema)) != "Ifc4x3_add2") { + SKIP("fixture is authored for IFC4X3_ADD2"); + } + CHECK(count_geo42_converting(RAKED_SECTIONED_SOLID_SPF) == 0); +} + +TEST_CASE("IfcSectionedSolidHorizontal follows a non-axis-aligned directrix", "[ifcgeom][infra-sweep]") { + if (std::string(STRINGIFY(IfcSchema)) != "Ifc4x3_add2") { + SKIP("fixture is authored for IFC4X3_ADD2"); + } + const auto result = convert_sectioned_solid(AXIS_ALIGNED_SECTIONED_SOLID_SPF); + CHECK(result.geo42_count == 0); + REQUIRE(result.produced_brep); + // Plan projection is the top plus the bottom of the slab, 2 x width x length; + // a collapsed sweep (the pre-fix behaviour) is nowhere near this. + CHECK(result.projected_area_z == Catch::Approx(2.0 * 12.0 * 60.0).margin(2.0)); +} + TEST_CASE("IfcGeom C++ fixture creates walls below and above the void limit", "[ifcgeom][voids]") { hierarchy_helper below_limit_file; const auto below_limit_wall = create_wall_with_voids(below_limit_file, MAX_VOIDS - 1); diff --git a/src/ifcopenshell-python/test/test_create_shape.py b/src/ifcopenshell-python/test/test_create_shape.py index 4118685dde..f9f9fd1af5 100644 --- a/src/ifcopenshell-python/test/test_create_shape.py +++ b/src/ifcopenshell-python/test/test_create_shape.py @@ -1,10 +1,12 @@ import functools import itertools +import math import multiprocessing import operator import os from typing import get_args +import numpy as np import pytest import ifcopenshell @@ -270,6 +272,276 @@ def test_logging(): ] +class TestSectionedSolidHorizontalRakedEndCut(test.bootstrap.IFC4X3): + """An IfcSectionedSolidHorizontal whose two IfcAxis2PlacementLinear cross + section positions use direction vectors inconsistently -- one raked + RefDirection + width scale, the other plain -- must still loft a uniform + prism (raked at one end, square at the other), not a wedge, and must not + log GEO 42.""" + + def _build_wingwall(self, theta): + f = self.file + ifcopenshell.api.root.create_entity(f, ifc_class="IfcProject", name="Test") + ctx = ifcopenshell.api.context.add_context(f, context_type="Model") + body = ifcopenshell.api.context.add_context( + f, context_type="Model", context_identifier="Body", target_view="MODEL_VIEW", parent=ctx + ) + + self.L, self.w, self.h = 40.0, 20.0 / 12.0, 6.0 + width_scale = 1.0 / math.cos(theta) + + directrix = f.createIfcPolyline( + Points=[f.createIfcCartesianPoint((0.0, 0.0, 0.0)), f.createIfcCartesianPoint((self.L, 0.0, 0.0))] + ) + + def rect(half_width): + coords = ( + (-half_width, 0.0), + (half_width, 0.0), + (half_width, self.h), + (-half_width, self.h), + (-half_width, 0.0), + ) + return f.createIfcArbitraryClosedProfileDef( + ProfileType="AREA", + OuterCurve=f.createIfcIndexedPolyCurve( + Points=f.createIfcCartesianPointList2D(coords), Segments=None, SelfIntersect=False + ), + ) + + axis = f.createIfcDirection((0.0, 0.0, 1.0)) + # Rakes the near end cap about the vertical axis by theta while keeping + # the section's width axis (Axis x RefDirection) perpendicular extent at w. + ref_direction = f.createIfcDirection((math.cos(theta), -math.sin(theta), 0.0)) + + def location(distance_along): + return f.createIfcPointByDistanceExpression( + DistanceAlong=f.createIfcLengthMeasure(distance_along), BasisCurve=directrix + ) + + near = f.createIfcAxis2PlacementLinear(Location=location(0.0), Axis=axis, RefDirection=ref_direction) + far = f.createIfcAxis2PlacementLinear(Location=location(self.L), Axis=axis) + + solid = f.createIfcSectionedSolidHorizontal( + Directrix=directrix, + CrossSections=[rect(0.5 * self.w * width_scale), rect(0.5 * self.w)], + CrossSectionPositions=[near, far], + ) + return f.createIfcShapeRepresentation( + ContextOfItems=body, + RepresentationIdentifier="Body", + RepresentationType="AdvancedSweptSolid", + Items=[solid], + ) + + @pytest.mark.skipif( + not ifcopenshell.geom.has_geometry_library("opencascade"), reason="requires the OpenCASCADE kernel" + ) + def test_raked_end_keeps_uniform_perpendicular_thickness(self): + theta = math.radians(25.0 + 22.0 / 60.0 + 58.0 / 3600.0) + representation = self._build_wingwall(theta) + + logger = ifcopenshell.logger() + logger.output_format(logger.FMT_INMEMORY) + settings = ifcopenshell.geom.settings() + settings.set("use-world-coords", True) + geometry = ifcopenshell.geom.create_shape(settings, representation, logger=logger) + + assert "GEO42" not in [msg.code for msg in logger] + + v = ifcopenshell.util.shape.get_vertices(geometry) + xs = v[:, 0] + + # Uniform thickness perpendicular to the (X aligned) directrix. The whole + # solid spans exactly w in Y; a wedge (the pre fix behaviour) does not. + assert np.ptp(v[:, 1]) == pytest.approx(self.w, abs=1e-4) + for x0 in np.linspace(3.0, self.L - 3.0, 12): + slab = v[np.abs(xs - x0) < 1.0] + if len(slab) < 4: + continue + assert np.ptp(slab[:, 1]) == pytest.approx(self.w, abs=1e-4), f"thickness at x={x0:.1f}" + + # Far end square to the directrix (all vertices at x == L), near end raked + # about the vertical axis by exactly theta (its extreme vertex sits at + # -w/2 * tan(theta) along the directrix). + assert xs.max() == pytest.approx(self.L, abs=1e-4) + assert np.ptp(v[xs > xs.max() - 1e-4][:, 0]) < 1e-4 + assert xs.min() == pytest.approx(-0.5 * self.w * math.tan(theta), abs=1e-4) + + +class TestSectionedSolidHorizontalHonoursAxis(test.bootstrap.IFC4X3): + """An IfcSectionedSolidHorizontal whose cross section placements carry an + explicit Axis = (0,0,1) but no RefDirection must sweep the profile with its + Y axis on that Axis and its normal on the directrix tangent (buildingSMART + IFC4.x-IF #147) -- i.e. following the curve. Regression test for the case + where the directrix does not run along the global X axis: the profile used + to be placed with a fixed [e_y | e_z | e_x] permutation that ignored the + directrix, collapsing the swept solid (a road pavement running north-south + would come out a sliver a few centimetres wide).""" + + def _build(self, theta=0.0): + f = self.file + ifcopenshell.api.root.create_entity(f, ifc_class="IfcProject", name="Test") + ctx = ifcopenshell.api.context.add_context(f, context_type="Model") + body = ifcopenshell.api.context.add_context( + f, context_type="Model", context_identifier="Body", target_view="MODEL_VIEW", parent=ctx + ) + + self.L, self.width, self.height = 60.0, 12.0, 0.5 + + # Directrix runs along +Y, not +X. + directrix = f.createIfcPolyline( + Points=[f.createIfcCartesianPoint((0.0, 0.0, 0.0)), f.createIfcCartesianPoint((0.0, self.L, 0.0))] + ) + + def rect(half_width): + coords = ( + (-half_width, 0.0), + (half_width, 0.0), + (half_width, self.height), + (-half_width, self.height), + (-half_width, 0.0), + ) + return f.createIfcArbitraryClosedProfileDef( + ProfileType="AREA", + OuterCurve=f.createIfcIndexedPolyCurve( + Points=f.createIfcCartesianPointList2D(coords), Segments=None, SelfIntersect=False + ), + ) + + axis = f.createIfcDirection((0.0, 0.0, 1.0)) + + def location(distance_along): + return f.createIfcPointByDistanceExpression( + DistanceAlong=f.createIfcLengthMeasure(distance_along), BasisCurve=directrix + ) + + near = f.createIfcAxis2PlacementLinear(Location=location(0.0), Axis=axis) + if theta: + # Directrix tangent is +Y; rake the far cap about the vertical Axis. + ref_direction = f.createIfcDirection((math.sin(theta), math.cos(theta), 0.0)) + far = f.createIfcAxis2PlacementLinear(Location=location(self.L), Axis=axis, RefDirection=ref_direction) + far_half = 0.5 * self.width / math.cos(theta) + else: + far = f.createIfcAxis2PlacementLinear(Location=location(self.L), Axis=axis) + far_half = 0.5 * self.width + + solid = f.createIfcSectionedSolidHorizontal( + Directrix=directrix, + CrossSections=[rect(0.5 * self.width), rect(far_half)], + CrossSectionPositions=[near, far], + ) + return f.createIfcShapeRepresentation( + ContextOfItems=body, + RepresentationIdentifier="Body", + RepresentationType="AdvancedSweptSolid", + Items=[solid], + ) + + def _shape(self, representation): + logger = ifcopenshell.logger() + logger.output_format(logger.FMT_INMEMORY) + settings = ifcopenshell.geom.settings() + settings.set("use-world-coords", True) + geometry = ifcopenshell.geom.create_shape(settings, representation, logger=logger) + assert "GEO42" not in [msg.code for msg in logger] + return ifcopenshell.util.shape.get_vertices(geometry) + + @pytest.mark.skipif( + not ifcopenshell.geom.has_geometry_library("opencascade"), reason="requires the OpenCASCADE kernel" + ) + def test_square_run_follows_the_directrix(self): + v = self._shape(self._build()) + # width across the road -> world X, crown -> world Z, length -> world Y. + assert np.ptp(v[:, 0]) == pytest.approx(self.width, abs=1e-4) + assert np.ptp(v[:, 1]) == pytest.approx(self.L, abs=1e-4) + assert np.ptp(v[:, 2]) == pytest.approx(self.height, abs=1e-4) + + @pytest.mark.skipif( + not ifcopenshell.geom.has_geometry_library("opencascade"), reason="requires the OpenCASCADE kernel" + ) + def test_raked_far_end_on_a_non_axis_aligned_directrix(self): + theta = math.radians(25.0 + 22.0 / 60.0 + 58.0 / 3600.0) + v = self._shape(self._build(theta)) + ys = v[:, 1] + # Uniform width perpendicular to the directrix, square near end, raked far end. + assert np.ptp(v[:, 0]) == pytest.approx(self.width, abs=1e-4) + for y0 in np.linspace(5.0, self.L - 5.0, 10): + slab = v[np.abs(ys - y0) < 1.0] + if len(slab) >= 4: + assert np.ptp(slab[:, 0]) == pytest.approx(self.width, abs=1e-4), f"width at y={y0:.1f}" + assert ys.min() == pytest.approx(0.0, abs=1e-4) + assert np.ptp(v[ys < ys.min() + 1e-4][:, 1]) < 1e-4 + assert ys.max() == pytest.approx(self.L + 0.5 * self.width * math.tan(theta), abs=1e-4) + + +class TestSectionedSolidHorizontalOffsetUnits(test.bootstrap.IFC4X3): + """IfcPointByDistanceExpression.OffsetLateral / OffsetVertical are + IfcLengthMeasure and must be scaled by the model length unit, exactly like + DistanceAlong. Regression test: in a millimetre model a 3000 mm lateral + offset must move the swept solid 3 m, not 3000 m.""" + + @pytest.mark.skipif( + not ifcopenshell.geom.has_geometry_library("opencascade"), reason="requires the OpenCASCADE kernel" + ) + def test_lateral_offset_is_scaled_by_the_length_unit(self): + f = self.file + ifcopenshell.api.root.create_entity(f, ifc_class="IfcProject", name="Test") + unit = ifcopenshell.api.unit.add_si_unit(f, unit_type="LENGTHUNIT", prefix="MILLI") + ifcopenshell.api.unit.assign_unit(f, units=[unit]) + ctx = ifcopenshell.api.context.add_context(f, context_type="Model") + body = ifcopenshell.api.context.add_context( + f, context_type="Model", context_identifier="Body", target_view="MODEL_VIEW", parent=ctx + ) + + length_mm, width_mm, offset_mm = 20000.0, 4000.0, 3000.0 + directrix = f.createIfcPolyline( + Points=[f.createIfcCartesianPoint((0.0, 0.0, 0.0)), f.createIfcCartesianPoint((length_mm, 0.0, 0.0))] + ) + hw = 0.5 * width_mm + profile = f.createIfcArbitraryClosedProfileDef( + ProfileType="AREA", + OuterCurve=f.createIfcIndexedPolyCurve( + Points=f.createIfcCartesianPointList2D(((-hw, 0.0), (hw, 0.0), (hw, 500.0), (-hw, 500.0), (-hw, 0.0))), + Segments=None, + SelfIntersect=False, + ), + ) + axis = f.createIfcDirection((0.0, 0.0, 1.0)) + + def position(distance_along): + return f.createIfcAxis2PlacementLinear( + Location=f.createIfcPointByDistanceExpression( + DistanceAlong=f.createIfcLengthMeasure(distance_along), + OffsetLateral=offset_mm, + BasisCurve=directrix, + ), + Axis=axis, + ) + + solid = f.createIfcSectionedSolidHorizontal( + Directrix=directrix, + CrossSections=[profile, profile], + CrossSectionPositions=[position(0.0), position(length_mm)], + ) + representation = f.createIfcShapeRepresentation( + ContextOfItems=body, + RepresentationIdentifier="Body", + RepresentationType="AdvancedSweptSolid", + Items=[solid], + ) + + settings = ifcopenshell.geom.settings() + settings.set("use-world-coords", True) + geometry = ifcopenshell.geom.create_shape(settings, representation) + v = ifcopenshell.util.shape.get_vertices(geometry) # metres + + # Directrix +X, Axis +Z -> profile local x is world +Y; +OffsetLateral shifts there. + assert np.ptp(v[:, 0]) == pytest.approx(length_mm / 1000.0, abs=1e-4) + assert np.ptp(v[:, 1]) == pytest.approx(width_mm / 1000.0, abs=1e-4) + assert v[:, 1].mean() == pytest.approx(offset_mm / 1000.0, abs=1e-3) + + if __name__ == "__main__": import pytest