From 6510291f219611e7ae88a2348388ab267e84fbbf Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Thu, 21 Sep 2023 15:48:12 -0700 Subject: [PATCH 1/7] Updated build scripts to use boost_1.79.0 --- win/build-deps.cmd | 2 +- win/run-cmake.bat | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/win/build-deps.cmd b/win/build-deps.cmd index 07acffd6ec..9d5ac6f07a 100644 --- a/win/build-deps.cmd +++ b/win/build-deps.cmd @@ -114,7 +114,7 @@ if "%CMAKE_VERSION%" LSS "cmake version 3.11.4" ( ) :: NOTE Boost < 1.64 doesn't work without tricks if the user has only VS 2017 installed and no earlier versions. -set BOOST_VERSION=1.74.0 +set BOOST_VERSION=1.79.0 :: Version string with underscores instead of dots. set BOOST_VER=%BOOST_VERSION:.=_% diff --git a/win/run-cmake.bat b/win/run-cmake.bat index 92741419bc..b8c81f03f9 100755 --- a/win/run-cmake.bat +++ b/win/run-cmake.bat @@ -80,7 +80,7 @@ IF NOT EXIST ..\%BUILD_DIR%. mkdir ..\%BUILD_DIR% pushd ..\%BUILD_DIR% :: tfk: todo remove duplication -set BOOST_VERSION=1.74.0 +set BOOST_VERSION=1.79.0 set BOOST_VER=%BOOST_VERSION:.=_% set BOOST_ROOT=%DEPS_DIR%\boost_%BOOST_VER% From 60723cff2b6af1fadd464b44899b58fa85b68500 Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Thu, 21 Sep 2023 15:47:40 -0700 Subject: [PATCH 2/7] Fixed problem with call to boost::mpl::for_each --- src/ifcgeom/mapping/IfcCurveSegment.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/ifcgeom/mapping/IfcCurveSegment.cpp b/src/ifcgeom/mapping/IfcCurveSegment.cpp index 61fc913c15..9d13729d8d 100644 --- a/src/ifcgeom/mapping/IfcCurveSegment.cpp +++ b/src/ifcgeom/mapping/IfcCurveSegment.cpp @@ -121,7 +121,7 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcCurveSegment* inst) { static int NUM_SEGMENTS = 64; curve_segment_evaluator cse(length_unit_, inst->ParentCurve(), inst->SegmentStart(), inst->SegmentLength()); - boost::mpl::for_each>(cse); + boost::mpl::for_each>(std::ref(cse)); std::vector polygon; From 4d3a21de7703ba2fa8ae1eca571a16d99c6e0c29 Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Fri, 22 Sep 2023 13:49:18 -0700 Subject: [PATCH 3/7] Fixes cast to IfcCurveSegment for access to ParentCurve --- src/ifcgeom/mapping/IfcCompositeCurve.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/ifcgeom/mapping/IfcCompositeCurve.cpp b/src/ifcgeom/mapping/IfcCompositeCurve.cpp index 7ed02dfabb..136953322a 100644 --- a/src/ifcgeom/mapping/IfcCompositeCurve.cpp +++ b/src/ifcgeom/mapping/IfcCompositeCurve.cpp @@ -66,7 +66,7 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcCompositeCurve* inst) { } #ifdef SCHEMA_HAS_IfcCurveSegment else if (segment->as()) { - auto crv = map(segment->as()->ParentCurve()); + auto crv = map(segment->as()->ParentCurve()); for (auto& s : taxonomy::cast(crv)->children) { loop->children.push_back(s); } From 9ac99b6f9775d17d5c94660230cd738605936497 Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Fri, 22 Sep 2023 13:52:48 -0700 Subject: [PATCH 4/7] Implements mapping of IfcCurveSegment to taxonomy::geometry --- src/ifcgeom/mapping/IfcCompositeCurve.cpp | 2 +- src/ifcgeom/mapping/IfcCurveSegment.cpp | 209 ++++++++++++++++++++-- 2 files changed, 198 insertions(+), 13 deletions(-) diff --git a/src/ifcgeom/mapping/IfcCompositeCurve.cpp b/src/ifcgeom/mapping/IfcCompositeCurve.cpp index 136953322a..1124cd083e 100644 --- a/src/ifcgeom/mapping/IfcCompositeCurve.cpp +++ b/src/ifcgeom/mapping/IfcCompositeCurve.cpp @@ -66,7 +66,7 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcCompositeCurve* inst) { } #ifdef SCHEMA_HAS_IfcCurveSegment else if (segment->as()) { - auto crv = map(segment->as()->ParentCurve()); + auto crv = map(segment->as()); for (auto& s : taxonomy::cast(crv)->children) { loop->children.push_back(s); } diff --git a/src/ifcgeom/mapping/IfcCurveSegment.cpp b/src/ifcgeom/mapping/IfcCurveSegment.cpp index 9d13729d8d..1de6bfaf5f 100644 --- a/src/ifcgeom/mapping/IfcCurveSegment.cpp +++ b/src/ifcgeom/mapping/IfcCurveSegment.cpp @@ -33,6 +33,8 @@ typedef boost::mpl::vector< #ifdef SCHEMA_HAS_IfcClothoid , IfcSchema::IfcClothoid #endif + , IfcSchema::IfcPolyline + , IfcSchema::IfcCircle > curve_seg_types; class curve_segment_evaluator { @@ -65,32 +67,189 @@ public: // Then initialize Function(double) -> Vector3, by means of IfcCurve subtypes void operator()(IfcSchema::IfcClothoid* c) { // @todo verify - auto L = start_ + length_; + auto sign = [](double v)->int{return v < 0 ? -1 : (0 < v ? 1 : 0); }; + auto sign_s = sign(start_); + auto sign_l = sign(length_); + double L = 0; + if (sign_s == 0) L = fabs(length_); + else if (sign_s == sign_l) L = fabs(start_ + length_); + else L = fabs(start_); + auto A = c->ClothoidConstant(); auto R = A * A / L; - auto RL = R * L; + auto RL = (A < 0 ? -1.0 : 1.0) * R * L; - eval_ = [RL](double u) { + auto position = c->Position(); + auto placement = position->as(); + auto ref_direction = placement->RefDirection(); + double theta = 0.0; // angle the circle's placement X-axis makes with respect to global X axis + if (ref_direction) + { + auto dr = ref_direction->DirectionRatios(); + auto dx = dr[0]; + auto dy = dr[1]; + theta = atan2(dy, dx); + } + + auto C = placement->Location(); + if (!C->as()) + { + throw std::runtime_error("Only IfcCartesianPoint is supported for center of IfcCircle"); + // @todo add support for other IfcPoint subtypes + } + auto Cx = C->as()->Coordinates()[0]; + auto Cy = C->as()->Coordinates()[1]; + + eval_ = [RL,Cx,Cy,theta](double u) { + // coordinate along clothoid is local coordinates auto xterm_1 = u; auto xterm_2 = std::pow(u, 5) / (40 * std::pow(RL, 2)); auto xterm_3 = std::pow(u, 9) / (3456 * std::pow(RL, 4)); auto xterm_4 = std::pow(u, 13) / (599040 * std::pow(RL, 6)); - auto x = xterm_1 - xterm_2 + xterm_3 - xterm_4; + auto xl = xterm_1 - xterm_2 + xterm_3 - xterm_4; auto yterm_1 = std::pow(u, 3) / (6 * RL); auto yterm_2 = std::pow(u, 7) / (336 * std::pow(RL, 3)); auto yterm_3 = std::pow(u, 11) / (42240 * std::pow(RL, 5)); auto yterm_4 = std::pow(u, 15) / (9676800 * std::pow(RL, 7)); - auto y = yterm_1 - yterm_2 + yterm_3 - yterm_4; + auto yl = yterm_1 - yterm_2 + yterm_3 - yterm_4; - return Eigen::Vector3d(x, y, 0.); + // transform point into clothoid's coodinate system + auto x = xl * cos(theta) - yl * sin(theta) + Cx; + auto y = xl * sin(theta) + yl * cos(theta) + Cy; + return Eigen::Vector3d(x, y, 0.0); }; } #endif - // Another IfcCurve subtype - void operator()(IfcSchema::IfcLine*) { - throw std::runtime_error("not implemented"); + void operator()(IfcSchema::IfcCircle* c) + { + auto R = c->Radius(); + + auto position = c->Position(); + auto placement = position->as(); + auto ref_direction = placement->RefDirection(); + double theta = 0.0; // angle the circle's placement X-axis makes with respect to global X axis + if (ref_direction) + { + auto dr = ref_direction->DirectionRatios(); + auto dx = dr[0]; + auto dy = dr[1]; + theta = atan2(dy, dx); + } + + // center of circle location + auto C = placement->Location(); + if (!C->as()) + { + throw std::runtime_error("Only IfcCartesianPoint is supported for center of IfcCircle"); + // @todo add support for other IfcPoint subtypes + } + auto Cx = C->as()->Coordinates()[0]; + auto Cy = C->as()->Coordinates()[1]; + + eval_ = [R, Cx, Cy, theta](double u) + { + auto angle = u / R; // angle subtended by arc length u + + // compute point on circle centered at (0,0) with x-axis horizontal and y-axis vertical + auto xl = R * cos(angle); + auto yl = R * sin(angle); + + // transform point into circle's coodinate system + auto x = xl * cos(theta) - yl * sin(theta) + Cx; + auto y = xl * sin(theta) + yl * cos(theta) + Cy; + return Eigen::Vector3d(x, y, 0.0); + }; + } + + void operator()(IfcSchema::IfcPolyline* pl) + { + struct Range + { + double u_start; + double u_end; + std::function compare; + bool operator<(const Range& r) const { return u_start < r.u_start; } + }; + using Function = std::function(double u)>; + std::map fns; + + auto p = pl->Points(); + if (p->size() < 2) + { + throw std::runtime_error("invalid polyline - must have at least 2 points"); // this should never happen, but just in case it does + } + + auto std_compare = [](double u_start, double u, double u_end) {return u_start <= u && u < u_end; }; + auto end_compare = [](double u_start, double u, double u_end) {return u_start <= u && u <= (u_end+0.001); }; + + auto iter = p->begin(); + auto end = p->end(); + auto last = std::prev(end); + auto p1 = *(iter++); + auto u = 0.0; + for (; iter != end; iter++) + { + auto p2 = *iter; + + auto p1x = p1->Coordinates()[0]; + auto p1y = p1->Coordinates()[1]; + + auto p2x = p2->Coordinates()[0]; + auto p2y = p2->Coordinates()[1]; + + auto dx = p2x - p1x; + auto dy = p2y - p1y; + auto l = sqrt(dx * dx + dy * dy); + if (l == 0.0) + { + // @todo use closeness tolerance instead of absolute 0.0 + throw std::runtime_error("invalid polyline - points must not be coincident"); + } + + dx /= l; + dy /= l; + + auto fn = [p1x, p1y, dx, dy](double u) { return std::make_pair(p1x + u * dx, p1y + u * dy); }; + + fns.insert(std::make_pair(Range{ u, u + l,iter == last ? end_compare : std_compare }, fn)); + + p1 = p2; + u = u + l; + } + + eval_ = [fns](double u) { + auto iter = std::find_if(fns.cbegin(), fns.cend(), [=](const auto& fn) + { + auto [u_start, u_end, compare] = fn.first; + return compare(u_start, u, u_end); + }); + + if (iter == fns.end()) throw std::runtime_error("invalid distance from start"); // this should never happen, but just in case it does + + auto [u_start, u_end, compare] = iter->first; + auto [x,y] = (iter->second)(u - u_start); // (u - u_start) is distance from start of this segment of the polyline + return Eigen::Vector3d(x, y, 0); + }; + } + + void operator()(IfcSchema::IfcLine* l) { + auto s = l->Pnt(); + auto c = s->Coordinates(); + auto v = l->Dir(); + auto dr = v->Orientation()->DirectionRatios(); + auto m = v->Magnitude(); + auto px = c[0]; + auto py = c[1]; + auto dx = dr[0] / m; + auto dy = dr[1] / m; + + eval_ = [px, py, dx, dy](double u) { + auto x = px + u * dx; + auto y = py + u * dy; + return Eigen::Vector3d(x, y, 0); + }; } // Take the boost::type value from mpl::for_each and test it against our curve instance @@ -118,6 +277,7 @@ public: taxonomy::ptr mapping::map_impl(const IfcSchema::IfcCurveSegment* inst) { // @todo fixed number of segments or fixed interval? // @todo placement + // @todo figure out what to do with the zero length segments at the end of compound curves static int NUM_SEGMENTS = 64; curve_segment_evaluator cse(length_unit_, inst->ParentCurve(), inst->SegmentStart(), inst->SegmentLength()); @@ -125,9 +285,34 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcCurveSegment* inst) { std::vector polygon; - for (int i = 0; i <= NUM_SEGMENTS; ++i) { - auto p = cse(cse.length() * i / NUM_SEGMENTS); - polygon.push_back(taxonomy::make(p(0), p(1), p(2))); + auto placement = inst->Placement(); + auto location = placement->Location(); + auto Cx = location->as()->Coordinates()[0]; + auto Cy = location->as()->Coordinates()[1]; + auto ref_dir = placement->as()->RefDirection(); + auto dx = ref_dir->DirectionRatios()[0]; + auto dy = ref_dir->DirectionRatios()[1]; + auto angle = atan2(dy, dx); + + auto cos_angle = cos(angle); + auto sin_angle = sin(angle); + + auto length = cse.length(); + if (0.001 < fabs(length)) + { + for (int i = 0; i <= NUM_SEGMENTS; ++i) { + auto u = length * i / NUM_SEGMENTS; + + auto p = cse(u); + auto xl = p(0); + auto yl = p(1); + auto z = p(2); + + auto x = xl * cos_angle - yl * sin_angle + Cx; + auto y = xl * sin_angle + yl * cos_angle + Cy; + + polygon.push_back(taxonomy::make(x,y,z)); + } } return polygon_from_points(polygon); From 2d14ad1667add4ac2a7e294d8b2a756bfc2a0822 Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Thu, 28 Sep 2023 11:51:26 -0700 Subject: [PATCH 5/7] Revises implementation of operator(IfcSchema::IfcClothoid*) to use trapezoid rule numeric integration. This is a more general approach than the Taylor Series approximation and is applicable to all the other spiral types. I stubbed out a functor for IfcSecondOrderPolynomialSpiral, but haven't tested it. The purpose is to show the concept. --- src/ifcgeom/mapping/IfcCurveSegment.cpp | 191 ++++++++++++++++++++---- 1 file changed, 161 insertions(+), 30 deletions(-) diff --git a/src/ifcgeom/mapping/IfcCurveSegment.cpp b/src/ifcgeom/mapping/IfcCurveSegment.cpp index 1de6bfaf5f..cd92d229bb 100644 --- a/src/ifcgeom/mapping/IfcCurveSegment.cpp +++ b/src/ifcgeom/mapping/IfcCurveSegment.cpp @@ -28,10 +28,38 @@ using namespace ifcopenshell::geometry; #include #include +// @todo use std::numbers::pi when upgrading to C++ 20 +#define PI 3.1415926535897932384626433832795 + + +namespace +{ + // trapezoid rule integration + // @todo is there a well established math library we can use instead of + // creating our own integrator? + double integrate(double a, double b, unsigned n, std::function fn) + { + double area = 0; + double h = (b - a) / (n + 1); + for (auto i = 0; i < n; i++) + { + auto x1 = a + h * i; + auto x2 = a + h * (i + 1); + auto f1 = fn(x1); + auto f2 = fn(x2); + area += h * (f1 + f2) / 2.0; + } + return area; + } +} + typedef boost::mpl::vector< IfcSchema::IfcLine #ifdef SCHEMA_HAS_IfcClothoid , IfcSchema::IfcClothoid +#endif +#if defined SCHEMA_HAS_IfcSecondOrderPolynomialSpiral + , IfcSchema::IfcSecondOrderPolynomialSpiral #endif , IfcSchema::IfcPolyline , IfcSchema::IfcCircle @@ -63,24 +91,80 @@ public: length_ = *le->as() * length_unit; } -#ifdef SCHEMA_HAS_IfcClothoid - // Then initialize Function(double) -> Vector3, by means of IfcCurve subtypes - void operator()(IfcSchema::IfcClothoid* c) { - // @todo verify - auto sign = [](double v)->int{return v < 0 ? -1 : (0 < v ? 1 : 0); }; - auto sign_s = sign(start_); - auto sign_l = sign(length_); +// Clothoid using Taylor Series approximation +//#ifdef SCHEMA_HAS_IfcClothoid +// // Then initialize Function(double) -> Vector3, by means of IfcCurve subtypes +// void operator()(IfcSchema::IfcClothoid* c) { +// // @todo verify +// auto sign = [](double v)->int{return v < 0 ? -1 : (0 < v ? 1 : 0); }; +// auto sign_s = sign(start_); +// auto sign_l = sign(length_); +// double L = 0; +// if (sign_s == 0) L = fabs(length_); +// else if (sign_s == sign_l) L = fabs(start_ + length_); +// else L = fabs(start_); +// +// auto A = c->ClothoidConstant(); +// auto R = A * A / L; +// auto RL = (A < 0 ? -1.0 : 1.0) * R * L; +// +// auto position = c->Position(); +// auto placement = position->as(); +// auto ref_direction = placement->RefDirection(); +// double theta = 0.0; // angle the circle's placement X-axis makes with respect to global X axis +// if (ref_direction) +// { +// auto dr = ref_direction->DirectionRatios(); +// auto dx = dr[0]; +// auto dy = dr[1]; +// theta = atan2(dy, dx); +// } +// +// auto C = placement->Location(); +// if (!C->as()) +// { +// throw std::runtime_error("Only IfcCartesianPoint is supported for center of IfcCircle"); +// // @todo add support for other IfcPoint subtypes +// } +// auto Cx = C->as()->Coordinates()[0]; +// auto Cy = C->as()->Coordinates()[1]; +// +// eval_ = [RL,Cx,Cy,theta](double u) { +// // coordinate along clothoid is local coordinates +// auto xterm_1 = u; +// auto xterm_2 = std::pow(u, 5) / (40 * std::pow(RL, 2)); +// auto xterm_3 = std::pow(u, 9) / (3456 * std::pow(RL, 4)); +// auto xterm_4 = std::pow(u, 13) / (599040 * std::pow(RL, 6)); +// auto xl = xterm_1 - xterm_2 + xterm_3 - xterm_4; +// +// auto yterm_1 = std::pow(u, 3) / (6 * RL); +// auto yterm_2 = std::pow(u, 7) / (336 * std::pow(RL, 3)); +// auto yterm_3 = std::pow(u, 11) / (42240 * std::pow(RL, 5)); +// auto yterm_4 = std::pow(u, 15) / (9676800 * std::pow(RL, 7)); +// auto yl = yterm_1 - yterm_2 + yterm_3 - yterm_4; +// +// // transform point into clothoid's coodinate system +// auto x = xl * cos(theta) - yl * sin(theta) + Cx; +// auto y = xl * sin(theta) + yl * cos(theta) + Cy; +// return Eigen::Vector3d(x, y, 0.0); +// }; +// } +//#endif + + void set_spiral_functor(IfcSchema::IfcSpiral* s, std::function signX,std::function fnX, std::function signY, std::function fnY) + { + // determine the length of the spiral from the local origin to the end point + auto binary_sign = [](double v)->int {return v < 0 ? -1 : (0 < v ? 1 : 0); }; // returns -1, 0, or 1 + auto sign_s = binary_sign(start_); + auto sign_l = binary_sign(length_); double L = 0; - if (sign_s == 0) L = fabs(length_); - else if (sign_s == sign_l) L = fabs(start_ + length_); - else L = fabs(start_); + if (sign_s == 0) L = fabs(length_); // start_ is at zero so length_ is the L + else if (sign_s == sign_l) L = fabs(start_ + length_); // start_ and length_ are additive + else L = fabs(start_); // start_ and length_ are in opposite directions so start_ is furthest from the origin - auto A = c->ClothoidConstant(); - auto R = A * A / L; - auto RL = (A < 0 ? -1.0 : 1.0) * R * L; - - auto position = c->Position(); - auto placement = position->as(); + auto position = s->Position(); + auto placement = position->as(); // @todo Update, this could be IfcAxis2Placement2D or IfcAxisPlacement3D + if (!placement) { throw std::runtime_error("Only IfcAxis2Placement2D is supported right now"); } auto ref_direction = placement->RefDirection(); double theta = 0.0; // angle the circle's placement X-axis makes with respect to global X axis if (ref_direction) @@ -94,31 +178,78 @@ public: auto C = placement->Location(); if (!C->as()) { - throw std::runtime_error("Only IfcCartesianPoint is supported for center of IfcCircle"); + throw std::runtime_error("Only IfcCartesianPoint is supported right now"); // @todo add support for other IfcPoint subtypes } auto Cx = C->as()->Coordinates()[0]; auto Cy = C->as()->Coordinates()[1]; - eval_ = [RL,Cx,Cy,theta](double u) { - // coordinate along clothoid is local coordinates - auto xterm_1 = u; - auto xterm_2 = std::pow(u, 5) / (40 * std::pow(RL, 2)); - auto xterm_3 = std::pow(u, 9) / (3456 * std::pow(RL, 4)); - auto xterm_4 = std::pow(u, 13) / (599040 * std::pow(RL, 6)); - auto xl = xterm_1 - xterm_2 + xterm_3 - xterm_4; + eval_ = [L, Cx, Cy, theta, signX, fnX, signY, fnY](double u) { + // integration limits, integrate from a to b + // from 8.9.3.19.1, integration limits are 0.0 to u where u is a normalized parameter + auto a = 0.0; + auto b = fabs(u / L); - auto yterm_1 = std::pow(u, 3) / (6 * RL); - auto yterm_2 = std::pow(u, 7) / (336 * std::pow(RL, 3)); - auto yterm_3 = std::pow(u, 11) / (42240 * std::pow(RL, 5)); - auto yterm_4 = std::pow(u, 15) / (9676800 * std::pow(RL, 7)); - auto yl = yterm_1 - yterm_2 + yterm_3 - yterm_4; + auto n = 10; // use 10 steps in the numeric integration + + auto xl = signX(u)*integrate(a, b, n, fnX); + auto yl = signY(u)*integrate(a, b, n, fnY); // transform point into clothoid's coodinate system auto x = xl * cos(theta) - yl * sin(theta) + Cx; auto y = xl * sin(theta) + yl * cos(theta) + Cy; return Eigen::Vector3d(x, y, 0.0); - }; + }; + } + +// Clothoid using numerical integration +#ifdef SCHEMA_HAS_IfcClothoid +// Then initialize Function(double) -> Vector3, by means of IfcCurve subtypes + void operator()(IfcSchema::IfcClothoid* c) { + + auto A = c->ClothoidConstant(); + + // the integration is for the +X, +Y quadrant - need to adjust the signs of the resulting X and Y values + // so that the results are in the correct quadrant. + // A > 0 and u > 0 -> +X, +Y + // A < 0 and u > 0 -> +X, -Y + // A > 0 and u < 0 -> -X, -Y + // A < 0 and u < 0 -> -X, +Y + // X depends only on u, Y depends on u and A. + auto sign = [](double v)->int {return v < 0 ? -1 : 1; }; // returns -1 or 1 + auto sign_x = [sign](double t) {return sign(t); }; + auto sign_y = [sign, A](double t) {return sign(t) == sign(A) ? 1.0 : -1.0; }; + auto fn_x = [A](double t)->double {return A * sqrt(PI) * cos(PI * A * t * t / (2 * fabs(A))); }; + auto fn_y = [A](double t)->double {return A * sqrt(PI) * sin(PI * A * t * t / (2 * fabs(A))); }; + + set_spiral_functor(c->as(), sign_x, fn_x, sign_y, fn_y); + } +#endif + +#ifdef SCHEMA_HAS_IfcSecondOrderPolynomialSpiral + void operator()(IfcSchema::IfcSecondOrderPolynomialSpiral* s) + { + // @todo verify - this is an example implementation of a different kind of spiral - lots of clean up needed + auto A0 = s->ConstantTerm(); + auto A1 = s->LinearTerm(); + auto A2 = s->QuadraticTerm(); + + auto theta = [A0, A1, A2](double t) + { + auto a0 = A0.has_value() ? t / A0.value() : 0.0; + auto a1 = A1.has_value() ? A1.value() * std::pow(t, 2) / (2 * fabs(std::pow(A1.value(), 3))) : 0.0; + auto a2 = std::pow(t, 3) / (3 * std::pow(A2, 3)); + return a0 + a1 + a2; + }; + + auto sign = [](double v)->int {return v < 0 ? -1 : 1; }; // returns -1 or 1 + auto sign_x = [sign](double t) {return sign(t); }; + auto sign_y = [sign](double t) {return sign(t); }; // @todo fix - not sure about sign_y yet, need to find some plots of this spiral + + auto fn_x = [theta](double t)->double {return cos(theta(t)); }; + auto fn_y = [theta](double t)->double {return sin(theta(t)); }; + + set_spiral_functor(s->as(), sign_x, fn_x, sign_y, fn_y); } #endif From fdc0c0cbd12f865c61a1869a6cd034d86f5a74f0 Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Thu, 28 Sep 2023 14:28:17 -0700 Subject: [PATCH 6/7] Minor revision to numerical integration code --- src/ifcgeom/mapping/IfcCurveSegment.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/ifcgeom/mapping/IfcCurveSegment.cpp b/src/ifcgeom/mapping/IfcCurveSegment.cpp index cd92d229bb..727032895d 100644 --- a/src/ifcgeom/mapping/IfcCurveSegment.cpp +++ b/src/ifcgeom/mapping/IfcCurveSegment.cpp @@ -40,11 +40,11 @@ namespace double integrate(double a, double b, unsigned n, std::function fn) { double area = 0; - double h = (b - a) / (n + 1); - for (auto i = 0; i < n; i++) + double h = (b - a) / n; + for (auto i = 1; i <= n; i++) { - auto x1 = a + h * i; - auto x2 = a + h * (i + 1); + auto x1 = a + h * (i - 1); + auto x2 = a + h * i; auto f1 = fn(x1); auto f2 = fn(x2); area += h * (f1 + f2) / 2.0; From d50dcadf5c18b98a3cd1c2036a460e97a427486b Mon Sep 17 00:00:00 2001 From: Richard Brice <37087370+RickBrice@users.noreply.github.com> Date: Thu, 28 Sep 2023 14:58:21 -0700 Subject: [PATCH 7/7] Implemented placement handling with transformation matrix. --- src/ifcgeom/mapping/IfcCurveSegment.cpp | 24 +++++------------------- 1 file changed, 5 insertions(+), 19 deletions(-) diff --git a/src/ifcgeom/mapping/IfcCurveSegment.cpp b/src/ifcgeom/mapping/IfcCurveSegment.cpp index 727032895d..7d81a29eed 100644 --- a/src/ifcgeom/mapping/IfcCurveSegment.cpp +++ b/src/ifcgeom/mapping/IfcCurveSegment.cpp @@ -416,17 +416,9 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcCurveSegment* inst) { std::vector polygon; - auto placement = inst->Placement(); - auto location = placement->Location(); - auto Cx = location->as()->Coordinates()[0]; - auto Cy = location->as()->Coordinates()[1]; - auto ref_dir = placement->as()->RefDirection(); - auto dx = ref_dir->DirectionRatios()[0]; - auto dy = ref_dir->DirectionRatios()[1]; - auto angle = atan2(dy, dx); - - auto cos_angle = cos(angle); - auto sin_angle = sin(angle); + // @todo - for some reason this isn't working, the matrix gets all messed up + //const auto& transformation_matrix = taxonomy::cast(map(inst->Placement()))->ccomponents(); + auto transformation_matrix = taxonomy::cast(map(inst->Placement()))->ccomponents(); auto length = cse.length(); if (0.001 < fabs(length)) @@ -435,14 +427,8 @@ taxonomy::ptr mapping::map_impl(const IfcSchema::IfcCurveSegment* inst) { auto u = length * i / NUM_SEGMENTS; auto p = cse(u); - auto xl = p(0); - auto yl = p(1); - auto z = p(2); - - auto x = xl * cos_angle - yl * sin_angle + Cx; - auto y = xl * sin_angle + yl * cos_angle + Cy; - - polygon.push_back(taxonomy::make(x,y,z)); + auto result = transformation_matrix * Eigen::Vector4d(p(0),p(1),p(2), 1.); + polygon.push_back(taxonomy::make(result(0),result(1),result(2))); } }