First attempt at taxonomy::implicit_item #3827

aothms · 2023-09-29T20:26:11Z

I've tried to create a prototype of this ``implicit_item'' idea I was describing. I've never seen such an amount of composition of lambdas in C++ code base. I'm not sure about it. Also, having to reverse engineer whether a curvesegment is part of a horizontal/vertical/cant and choose the appropriate return type deep down in this functor is probably less than ideal, but keeps the amount of code minimal. Curious to hear your feedback. I also took the liberty to replace integrate() with something from boost.math.

Added you as a member @RickBrice in order to assign you as a reviewer.

RickBrice

It seems the general approach will work. I've made some comments and suggestions for your consideration.

Not a big fan of using so many lambda functions either. They make the code hard to unit test.

RickBrice · 2023-09-30T02:43:19Z

src/ifcgeom/mapping/IfcCurveSegment.cpp

 			// 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 n = 10; // use 10 steps in the numeric integration
+			// @todo where to plug this in?


I think this is the 4th argument to trapezoidal

I don't think so. the docs say:

Since the routine is adaptive, step sizes are halved continuously until a tolerance is reached. In order to control this tolerance, simply call the routine with an additional argument

https://live.boost.org/doc/libs/1_68_0/libs/math/doc/html/math_toolkit/trapezoidal.html

So it's not just a fixed increment, but some limit to the number of refinement iterations. So conceptually it'd be more like 2^10 (big question mark). That's why I was a bit afraid to plug it in, but the default is 12. So I think we can just keep the defaults.

RickBrice · 2023-09-30T02:52:04Z

src/ifcgeom/mapping/IfcGradientCurve.cpp

+		}
+	}
+
+	// @todo does this really make sense?


Yes, this is typically how it's handled when there is a horizontal and vertical definition.

This doesn't work for CT 4.1.7.1.1.4 Alignment Geometry- Segments because it seems to require that the ParentCurve be 3D. There isn't independent horizontal and vertical components for this CT.

We are making another assumption here - sort of - maybe it's just a distinction we need to be mindful of.

For CT Alignment Geometry Segment, the representation is on the IfcCurveSegment whereas the other CTs put the representation on the IfcAlignment.

because it seems to require that the ParentCurve be 3D

Where do you see this constraint encoded? I think that would be an error. I briefly looked and IfcShapeRepresentation.RepresentationType=Segment seems to just select IfcCurveSegment without any constraints on the parent curve.

IFC is a community where clear choices are seldomly made. So I think the idea is you do everything. Individual representations for the segments is the most likely thing missing from these aspects. But you're right in the sense that due to the IFC schema flexibility, you never really know.

I've updated my diagram.

dot code

digraph { node [shape=rect] IfcAlignment -> IfcRelNests_0[dir=back] IfcRelNests_0 -> { IfcAlignmentHorizontal, IfcAlignmentVertical, IfcAlignmentCant } IfcAlignmentHorizontal -> IfcRelNests_1[dir=back] IfcRelNests_1 -> { IfcAlignmentSegment_00, IfcAlignmentSegment_01, IfcAlignmentSegment_02 } [minlen=5] IfcAlignmentVertical -> IfcRelNests_2[dir=back] IfcRelNests_2 -> { IfcAlignmentSegment_10, IfcAlignmentSegment_11, IfcAlignmentSegment_12 } [minlen=3] IfcAlignmentCant -> IfcRelNests_3[dir=back] IfcRelNests_3 -> { IfcAlignmentSegment_20, IfcAlignmentSegment_21, IfcAlignmentSegment_22 } [minlen=1] IfcAlignmentSegment_00 -> IfcAlignmentHorizontalSegment_0 IfcAlignmentSegment_01 -> IfcAlignmentHorizontalSegment_1 IfcAlignmentSegment_02 -> IfcAlignmentHorizontalSegment_2 IfcAlignmentSegment_10 -> IfcAlignmentVerticalSegment_0 IfcAlignmentSegment_11 -> IfcAlignmentVerticalSegment_1 IfcAlignmentSegment_12 -> IfcAlignmentVerticalSegment_2 IfcAlignmentSegment_20 -> IfcAlignmentCantSegment_0 IfcAlignmentSegment_21 -> IfcAlignmentCantSegment_1 IfcAlignmentSegment_22 -> IfcAlignmentCantSegment_2 IfcAlignmentCantSegment_0 [label=IfcAlignmentCantSegment] IfcAlignmentCantSegment_1 [label=IfcAlignmentCantSegment] IfcAlignmentCantSegment_2 [label=IfcAlignmentCantSegment] IfcAlignmentHorizontalSegment_0 [label=IfcAlignmentHorizontalSegment] IfcAlignmentHorizontalSegment_1 [label=IfcAlignmentHorizontalSegment] IfcAlignmentHorizontalSegment_2 [label=IfcAlignmentHorizontalSegment] IfcAlignmentSegment_00 [label=IfcAlignmentSegment] IfcAlignmentSegment_01 [label=IfcAlignmentSegment] IfcAlignmentSegment_02 [label=IfcAlignmentSegment] IfcAlignmentSegment_10 [label=IfcAlignmentSegment] IfcAlignmentSegment_11 [label=IfcAlignmentSegment] IfcAlignmentSegment_12 [label=IfcAlignmentSegment] IfcAlignmentSegment_20 [label=IfcAlignmentSegment] IfcAlignmentSegment_21 [label=IfcAlignmentSegment] IfcAlignmentSegment_22 [label=IfcAlignmentSegment] IfcAlignmentVerticalSegment_0 [label=IfcAlignmentVerticalSegment] IfcAlignmentVerticalSegment_1 [label=IfcAlignmentVerticalSegment] IfcAlignmentVerticalSegment_2 [label=IfcAlignmentVerticalSegment] IfcRelNests_0 [label=N, shape=circle] IfcRelNests_1 [label=N, shape=circle] IfcRelNests_2 [label=N, shape=circle] IfcRelNests_3 [label=N, shape=circle] IfcAlignment -> IfcProductDefinitionShape IfcProductDefinitionShape -> IfcShapeRepresentation IfcShapeRepresentation -> IfcSegmentedReferenceCurve IfcAlignment -> spacing1 [style=invis] spacing1 -> IfcProductDefinitionShape [style=invis] IfcProductDefinitionShape -> spacing2 [style=invis] spacing2 -> IfcShapeRepresentation [style=invis] IfcShapeRepresentation -> spacing3 [style=invis] spacing3 -> IfcSegmentedReferenceCurve [style=invis] { IfcAlignment IfcProductDefinitionShape IfcShapeRepresentation IfcSegmentedReferenceCurve spacing1 spacing2 spacing3 rank=same } IfcSegmentedReferenceCurve -> { IfcCurveSegment_20, IfcCurveSegment_21, IfcCurveSegment_22 } [minlen=5] IfcSegmentedReferenceCurve -> IfcGradientCurve IfcGradientCurve -> { IfcCurveSegment_10, IfcCurveSegment_11, IfcCurveSegment_12 } [minlen=6] IfcGradientCurve -> IfcCompositeCurve IfcCompositeCurve -> { IfcCurveSegment_00, IfcCurveSegment_01, IfcCurveSegment_02 } [minlen=7] spacing1 [style=invis, width=2] spacing2 [style=invis, width=2] spacing3 [style=invis, width=2] edge [constraint=false] IfcAlignmentSegment_20 -> IfcCurveSegment_20 [style=dashed, label=Representation] IfcAlignmentSegment_21 -> IfcCurveSegment_21 [style=dashed, label=Representation] IfcAlignmentSegment_22 -> IfcCurveSegment_22 [style=dashed, label=Representation] IfcAlignmentSegment_10 -> IfcCurveSegment_10 [style=dashed, label=Representation] IfcAlignmentSegment_11 -> IfcCurveSegment_11 [style=dashed, label=Representation] IfcAlignmentSegment_12 -> IfcCurveSegment_12 [style=dashed, label=Representation] IfcAlignmentSegment_00 -> IfcCurveSegment_00 [style=dashed, label=Representation] IfcAlignmentSegment_01 -> IfcCurveSegment_01 [style=dashed, label=Representation] IfcAlignmentSegment_02 -> IfcCurveSegment_02 [style=dashed, label=Representation] IfcCurveSegment_20 [label=IfcCurveSegment] IfcCurveSegment_21 [label=IfcCurveSegment] IfcCurveSegment_22 [label=IfcCurveSegment] IfcCurveSegment_10 [label=IfcCurveSegment] IfcCurveSegment_11 [label=IfcCurveSegment] IfcCurveSegment_12 [label=IfcCurveSegment] IfcCurveSegment_00 [label=IfcCurveSegment] IfcCurveSegment_01 [label=IfcCurveSegment] IfcCurveSegment_02 [label=IfcCurveSegment] }

Where do you see this constraint encoded?

Perhaps I over stated my assumption. There are 4 CTs, the first 3 decompose the 3D space into two-2D representations that are combined to get a 3D result. I assumed the 4th CT, Segments, models the 3D space directly. The ParentCurve could be in a 2D horizontal xy plane, but this would be the Horizontal CT. Not sure if it make sense that ParentCurve be in a vertical sz plane without a horizontal component.

I'm not sure what use case was imagined for the 4th CT but to me it looks like it would apply for as-built survey data which would be 3D points or segments

Ok, I interpreted this differently. In my view the segment geometry is not really meant to be interesting in itself, but more to provide a more granular coupling between the business logic and the geometry. So the concept template just expresses that and doesn't really relate to a usecase. For surveying usecases it is possible to express an IfcAlignment representation simply as as 3d polyline.

(Or another possible top-level representation item for alignment: https://ifc43-docs.standards.buildingsmart.org/IFC/RELEASE/IFC4x3/HTML/lexical/IfcOffsetCurveByDistances.htm)

RickBrice · 2023-09-30T02:59:43Z

src/ifcgeom/mapping/IfcGradientCurve.cpp

+#ifdef SCHEMA_HAS_IfcGradientCurve
+
+taxonomy::ptr mapping::map_impl(const IfcSchema::IfcGradientCurve* inst) {
+	auto horizontal = taxonomy::cast<taxonomy::piecewise_function>(map(inst->BaseCurve()));


Is there a way to get horizontal from a cache if it was previously mapped for a Footprint/Curve2D case or the case where multiple alignments use the same horizontal layout but different vertical profiles?

I will work on caching next. It's trivial, it's just a std::map<uint, taxonomy::item::ptr>. With again the caveat that we should be a bit careful considering immutability.

What if the map function returns a shared pointer to a constant taxonomy item? std::shared_ptr

The map function creates the taxonomy::item and all other places in the code are consumers of the taxonomy::item. Do consumers have a need to alter the item? Probably not.

This would enforce immutability.

Do consumers have a need to alter the item? Probably not.

There are one or two cases that still need to be worked out with orientation on things like IfcOrientedEdge, where the consumer does indeed alter item (minimally). (In OpenCascade they have their own shared pointer / handle, but orientation is basically stored on the pointer, not in the data.)

RickBrice · 2023-09-30T03:01:43Z

src/ifcgeom/mapping/IfcGradientCurve.cpp

+		return vec;
+	};
+
+	// @todo where do we get the startdistalong from @civilx64's code?


Do you get it from segment.start and segment.length?

https://github.com/civilx64/IfcOpenShell/blob/c263b4f50ebbdcdacc5d3f548c561a1fb24dcd13/src/ifcopenshell-python/ifcopenshell/ifcgeom/IfcAlignmentHorizontal.py#L206

When iterating segments I keep a private variable _length which is incremented by segment.length for each segment.

_length is short for "length along the alignment up to the beginning of this particular segment".

Sorry, I wasn't very clear, this was mostly a note to myself. This was about:

https://github.com/IfcOpenShell/IfcOpenShell/pull/3148/files#diff-4b6d2c91bd1d5ed55b20ad2d6023a46e93c7bdd011e3ecf4727c172e5ae436ecR325

On the business logic side we have StartDistAlong on vertical. But on the geometry side for the first vertical segment we have:

three parameters that can be used to independently change the start point of the vertical wrt to the horizontal:

The vertical curve start point (line Pnt in this case)

SegmentStart which is an offset over the vertical curve

Segment placement, which can be

linear placement over horizontal:

same horizontal (good, IfcPointByDistanceExpression.DistanceAlong becomes essentially StartDistAlong on vertical)

different horizontal (invalid?)

localplacement:

same as horizontal (can be ignored)

different from horizontal (in theory we could project the localplacement position onto the curve to find a length along the curve?)

I'm also really confused about the placement altogether. You can introduce gradient by means of IfcLine.Dir but also with the placement X-axis. Which dimensionality is a vertical alignment curve to begin with? Is it a curve in the XY-plane with that maps $X -> Y$ or $u -> Y$.

The attribute description of Placement is some really cryptic prose:

Placement in the context of the curve using this segment. As insertion point SegmentStart is the reference point of Placement. RefDirection of Placement also specifies the sense of the trimmed segment of ParentCurve. RefDirection is bound to the paramettrization sense of the segment.

RickBrice · 2023-09-30T03:07:57Z

src/ifcgeom/taxonomy.h

+
+	Eigen::VectorXd evaluate(double u) const {
+		// @todo optimize, assume monotonic evaluation and store last evaluated segment?
+		for (auto& s : spans) {


Structured binding makes this more readable

`for (auto& [start,fn] : spans)
{
if(u < start) return fn(u);

u -= start;

}`

Thanks. I'm very used to this in Python, but in C++ it's not really part of my muscle memory. Keep in mind it's $length, fn$ not start. But that's even better as this makes it more readable.

RickBrice · 2023-09-30T03:28:09Z

src/ifcgeom/taxonomy.h

+struct piecewise_function : public item {
+	DECLARE_PTR(piecewise_function)
+
+	std::vector<std::pair<double, std::function<Eigen::VectorXd(double u)>>> spans;


Maybe vector of tuple<double,function,function, function> where the 3 functions return xyz point, z-axis direction, and ref direction for distance along, u. For linear placement, if z-axis and ref direction are not explicitly provided, they are derived from the curve.

This is a prime example of overuse of lambda expressions. How can we do unit tests on the equations?

I see what you mean. That's also what I meant earlier with my comment. Hiding this composition in a lambda is felt less then ideal.

https://github.com/aothms/IfcOpenShell/blob/622a21cd41c4b93b5939be30adc325ad06c9e86a/src/ifcgeom/mapping/IfcGradientCurve.cpp#L49-L56

Instead of tuple<double,function,function, function>. How do you feel about:

// let's try to come up with better names enum curve_semantics { xy, z, cant }; std::map<curve_semantics, std::function<Eigen::VectorXd(double u)>>

then I would also propose not to merge this with piecewise_function but create another implicit_item subtype.

How do you feel about:

The 3 functions aren't xy, z, and cant. They are functions that compute the three items needed for linear placement. For distance along curve, linear placement needs xyz point, direction of z-axis and ref direction to determine x-axis.

We are addressing the xyz point only and your 'composite' idea handles the xy,z, and cant aspects needed to get xyz points.

The z-axis and ref direction are optional. When omitted (I think but can't seem to find it in the spec) z-axis is (0,0,1) and ref direction is tangent to the curve

Oh thanks for clarifying. In that case. Why don't we just return an Eigen::MatrixXd? which can potentially be just a column vector, doesn't need to be square. Or could be a 4x3 (or 4x4) representing rotation and translation.

Naively it seems more efficient to have a single function instead of three separate ones. I can imagine the code paths for these three functions would be rather similar. It's a matter of returning the 1st derivative in addition to the 0th basically (?).

I can imagine you are considering cases where you calculate data that is not needed. E.g for graphic display the rotation is not needed. What we can also do is pass the shape of the matrix we expect as an argument to the function. (or a bool to differentiate between two cases)

std::function<Eigen::MatrixXd(double u, const std::array<size_t>& dims)>

(I guess default arguments is not possible on a std::function)

I think an Eigen::MatrixXd would work fine. I only suggested 3-functions because that's the way my literal brain sees the solution, I don't have a lot of experience with the Eigen library, so it doesn't enter my thought process - I'll try to be better in that regard.

In general, there isn't a lot of overhead in computing extra parameters you don't need (z-axis and ref-direction) but this might not be the case as the number of linear placements scales up. Additionally, the z-axis and ref-direction only need to be computed if not provided for by the linear placement.

Ok. let's go with MatrixXd then. When we merge your PR #3826 I'll try to add it in.

RickBrice · 2023-09-30T13:47:54Z

src/ifcgeom/mapping/IfcCurveSegment.cpp

@@ -27,32 +27,11 @@ using namespace ifcopenshell::geometry;

 #include <boost/mpl/vector.hpp>
 #include <boost/mpl/for_each.hpp>
+#include <boost/math/quadrature/trapezoidal.hpp>


Boost math also has constants like pi

Changed. Thanks.

RickBrice · 2023-09-30T13:51:20Z

src/ifcgeom/mapping/IfcCurveSegment.cpp

+
+			eval_ = [px, py, dx, dy](double u) {
+				// @todo this can't be correct
+				auto z = u * dy;


auto z = py + u*dy;

RickBrice · 2023-09-30T13:58:00Z

src/ifcgeom/mapping/IfcCurveSegment.cpp

-			auto x = px + u * dx;
-			auto y = py + u * dy;
-			return Eigen::Vector3d(x, y, 0);
+		if (segment_type_ == ST_HORIZONTAL) {


The treatment of line assumes the line is either in the xy plane for horizontal or the sz plane for vertical. I don't think the ifc specification requires this.

CT 4.1.7.1.1.4 Alignment Geometry- Segments seems to require that the ParentCurve be 3D

sz plane = developed elevation view of the alignment profile. s = distance along curve, z = elevation.

Ah thanks. That's exactly what I was beginning to wonder earlier on. If we don't have any guidance on this in the spec it's really an omission that needs to be corrected. Luckily we also have the validation service where such constraints can be incrementally rolled out (cc @civilx64)

I'm not sure if such constraints are required but they would make implementation much easier and rational.

Consider a line in 3D space. A single IfcLine entity could be used for both the horizontal and vertical definitions. The projection of the line onto the xy and sz planes defines the horizontal and vertical geometry. This is essentially a diagonal line through opposite corners of a cube. That's fairly easy to deal with.

Now consider a clothoid in an arbitrary plane in 3D space. I really don't want to figure out its projections onto the xy and sz planes.

Additionally, I doubt the resulting geometry is actually used in practice since roadway and rail design predates computers - our predecessors engineers likely didn't make these calculations be hand.

RickBrice · 2023-09-30T14:07:19Z

src/ifcgeom/mapping/IfcCurveSegment.cpp

@@ -65,19 +44,25 @@ typedef boost::mpl::vector<
 	, IfcSchema::IfcCircle
 > curve_seg_types;

+enum segment_type_t {


Instead of this, coupled with the if-else if-else blocks for each IfcCurve type, would it be cleaner to have three different map functions for IfcCurveSegment? map_horizontal(), map_vertical(), map_cant(). These would get called from IfcCompositeCurve, IfcGradiantCurve, and IfcSegmentedReferenceCurve, respectively.

I would be okay with it as it, but admittedly at this point I'm just hanging on for dear life trying to keep up with the expert-level C++ exchanges between you & @aothms 🤣

Continuing to think about this, the original approach is better than my idea. My idea would have, for each curve type, three implementations in three different locations. @aothms Original idea has everything for each curve type together. That approach is more clear and maintainable.

I had the same doubts myself. I would have also preferred a more statically typed set of three tuples of types for example. I think the approach now is likely to result in the least lines of code, but it's indeed unclear what is implemented and supported and potentially risky in the sense that a parent curve can be used in an interpretation that's not supported. Combined with returning Eigen::VectorXd there is actually very little typing here to prevent issues. I would propose to keep as is while prototyping, but I think we might have to change it in the future.

civilx64

I agree with @RickBrice comments and have added feedback where I can. Hope it is useful.

aothms · 2023-10-01T09:53:30Z

Thanks for the comments both. I've added caching to this PR as well.

RickBrice · 2023-10-01T15:41:48Z

The framework is coming along nicely. One thing that seems to be missing is correcting numerical error using the start point of the next segment and dealing with the zero length last segment. See buildingSMART/IFC4.3.x-development#702.

I don't think providing a nextSegment parameter to the functor operator()(IfcCurveSegent) will get the job done since it is the mapped geometry that might need adjusting.

Also, the "start" of the next curve might actually be at the end. An example is an exit spiral to a horizontal curve. In the examples I was working with, the exit spiral starts at -Length and has a length of Length. This makes the spiral start at a point of finite radius (radius equal to that of the adjacent horizontal curve) and end at the point of infinite radius. In this case, the start point would have numerical error and the end point would be at (0,0) in the spiral coordinates. The start point of the next segment isn't very helpful in this case.

Any ideas on addressing this?

aothms · 2023-10-02T08:13:49Z

I don't think providing a nextSegment parameter to the functor operator()(IfcCurveSegent) will get the job done since it is the mapped geometry that might need adjusting.

For the STEP composite curve this is more or exactly what we do though:

A function that takes (Shape, Shape) as arguments

IfcOpenShell/src/ifcgeom_schema_agnostic/wire_builder.h

Line 56 in edaabf5

void operator()(const TopoDS_Shape& a, const TopoDS_Shape& b, bool last);

Loop over wire and apply to pairs of consecutive segments

IfcOpenShell/src/ifcgeom_schema_agnostic/wire_builder.h

Line 70 in edaabf5

fn(previous, current, false);

Compute distance and adjust end points or add intermediate edge if needed

IfcOpenShell/src/ifcgeom_schema_agnostic/wire_builder.cpp

Line 162 in edaabf5

double dist = p1.Distance(p2);

Including e.g refitting a circle

IfcOpenShell/src/ifcgeom_schema_agnostic/wire_builder.cpp

Line 66 in edaabf5

GC_MakeCircle mc(p1, p2, p3);

But that's also because in the BRep geometry in open cascade we have a bit different constraints.

I don't know what are the typical way to deal with this in infra geom implementations and what the expected deviation is. From earlier screenshots and discussion I seem to remember these deviations can be considerable, not like floating point rounding errors, but order of magnitude centimeters maybe.

First of all since we rely on polygon_from_points() we get a connected wire of points regardless of precision issues.

We could / should adapt polygon_from_points() to remove (near-) duplicate points within model tolerance. But that won't work for centimeter scale deviations (I guess).

Do you only need to correct the last point of the segment to eliminate the most obvious visual jump or do want to "interpolate this change over the entire segment", e.g:

segment_0 = [p0, p1, p2, p3]
segment_1 = [p4, p5, p6]
d = p4 - p3
lengths = [ ||q - p || for p, q in pairs(segment_0) ]
total_length = sum(lengths)
combined = [p0, p1 + (lengths[0] / total_length) * d, p2 + (lengths[1] / total_length) * d, p3 + d, p5, p6]

(or maybe not interpolate linearly, but quadratically, I'm really just making things up here).

The difficulty is indeed because of all our function composition we loose track of what actually constitutes a segment and we only have one domain. So we can change that. Or we can incorporate (any kind of this blending/interpolation) into the function composition

Eigen::VectorXd evaluate(double u) const {
    // @todo optimize, assume monotonic evaluation and store last evaluated segment?
    for (auto it = spans.begin(); it != spans.end(); +it) {
        if (u < it->first) {
            if (it + 1 == spans.end()) {
                return it->second(u);
            } else {
                double w = u / it->first;
                return (1. - w) * it->second(u) + w * (it + 1)->second(u - it->first);
            }
        }
        u -= it->first;
    }
}

I'm not too worried about the end segment. We just really need to be sure we call evaluate(length) for the piece_wise function, but since we no longer compose the full curve from individual polylines, but just a have single continuous domain I think this is handled more or less implicitly.

Also, the "start" of the next curve might actually be at the end. An example is an exit spiral to a horizontal curve. In the examples I was working with, the exit spiral starts at -Length and has a length of Length. This makes the spiral start at a point of finite radius (radius equal to that of the adjacent horizontal curve) and end at the point of infinite radius. In this case, the start point would have numerical error and the end point would be at (0,0) in the spiral coordinates. The start point of the next segment isn't very helpful in this case.

I see what you mean. So I think the interpolation approach could be adjusted to be applied on the other end as well. But what do you do if you have two connected clothoids? (I guess it's not incredibly common). Is it necessary then to store in our "spans" some sort of error estimate for the start and end of the span and use that to determine the interpolation. And then maybe if you interpolate clothoid to clothoid the weight at the end gets 0.5 so you average the begin and end point?

RickBrice · 2023-10-02T14:42:09Z

At this time, I'm not too concerned about how geometry corrections are applied. I think it is more important that there is a framework for making corrections and we can experiment with specific techniques based on actual geometric errors encountered as the work develops. Maybe the Strategy design pattern provides a good approach.

But what do you do if you have two connected clothoids?

Good question? Maybe split the difference like you suggested.

For bridges, which is where all my experience is, we rarely deal with sections of alignment that are more than a few hundred feet (most common bridges are highway crossings), longer bridges like river crossings tend to be straight. In this subset of alignment, there isn't much numeric error.

Perhaps @civilx64 can provide some insight for dealing with longer alignments.

RickBrice · 2023-10-02T22:25:18Z

src/ifcgeom/taxonomy.h

+				return fn(u);
+			}
+			u -= length;
+		}


not all paths return a value.

Hm yes indeed. I don't know what's reasonable here. I have also heard that sometimes vertical alignments actually extend beyond the horizontal in case of e.g linear horizontal segments and you'd be assumed to simply extrapolate the linear segment. Is that a thing? Should we accommodate for that? We could add a check whether it's the last segment and in case it's last don't check for u < length that would also save us some difficulties in dealing with the last segment, because due to rounding errors, evaluating the last 0-length segment will likely never happen in our current case. The only case not handled is negative u then. I don't have an opinion what to do (either throw an exception or extrapolate first segment beyond it's start).

I will need to dig into the spec again. My recollection is that the range of horizontal and vertical needs to be the same, but I want to confirm before proceeding with that assumption. In either case I will pursue a solution

Maybe also something for @civilx64 and the validation service. I can imagine at least a warning would be good to have in place in case of disagreement in range.

I reviewed the spec and didn't find anything specifically addressing this. I do remember seeing a use case or example specifically focusing on vertical range shorter than horizontal, which is common - and even recommended - in my professional experience.

I recommend the following for our implementation in IFCOS:

Vertical range can be shorter than horizontal at the beginning, end, or both limits of horizontal.

z values for horizontal (x, y) beyond the range of vertical shall be null

Vertical range beyond horizontal will be logged as a warning in IFCOS and the validation service.

@civilx64

z values for horizontal (x, y) beyond the range of vertical shall be null

How is this accomplished? The result of evaluating IfcGradientCurve or IfcSegmentedReferenceCurve is a tuple xyz of real numbers.

One option would be to use std::numeric_limits<double>::quiet_NaN(). It's maybe not very pretty, but should be supported everywhere since it's part of the IEEE spec.

RickBrice · 2023-10-02T22:58:17Z

src/ifcgeom/taxonomy.h

+	virtual item::ptr evaluate() const;
+
+	Eigen::VectorXd evaluate(double u) const {
+		// @todo optimize, assume monotonic evaluation and store last evaluated segment?


If span had start and end instead of length, this method could be optimized something like:

static std::tuple<double, double, std::function<Eigen::VectorXd(double u)>>* last_span = nullptr; if(last_span) { // there was a previous evaluated span auto& [start,end,fn] = *last_span; if(start <= u && u < end) return fn(u - start); // u is in the range for the previous span, evaluate it } // no last span, or u is not in the range of the last span.... find the span that covers u auto iter = std::find_if(spans.begin(),spans.end(),[u](auto& s){return u <= s.start && u < s.end;}); if(iter == spans.end()) { // deal with special case of u = s.end for last segment or error if u is not in the range of any spans} else { last_span = &(*iter); auto& [start,end,fn] = *iter; return fn(u - start); }

Sounds good except for one detail. static means shared by all instances of this class, also other alignment curves (so you might end up evaluating a span of another alignment). Instead I would add this as a member so that it's local to the alignment. You can mark it mutable so that it can be updated even in a function marked const. Would you add this to your PR #3826 then?

Yes, I can add this to my pr. Good point about static. I'll address it appropriately.

@aothms @civilx64 I'm have a little trouble getting this working. Since the concept is to be more efficient, I think it is more important to focus on getting the geometry we need to advance IFCOS. I can return to making this more efficient at a later date.

Absolutely agree. Same for things like hashing and all. I can work on that independently later on.

civilx64 · 2023-10-03T01:08:40Z

Maybe the Strategy design pattern provides a good approach.

I would agree.

As for numeric error, if I'm tracking correctly this happens on transition curves that are approximated by e.g. Taylor series expansion. My hunch is that this would be minimized by numeric integration with a sufficient number of steps / closure tolerance.

I've never had a reason to consider this professionally prior to this point. Typically in horizontal layout the designer's focus is on the tangent lines and arcs, with transitions "connecting the dots" so to speak.

aothms · 2023-10-05T09:44:03Z

Superseded by #3826 but let's keep it open until last comments are addressed

RickBrice · 2023-10-05T10:38:36Z

src/ifcgeom/mapping/mapping.i

@@ -113,6 +113,9 @@ BIND(IfcEdge);
 BIND(IfcEdgeLoop);
 BIND(IfcPolyline);
 BIND(IfcPolyLoop);
+#ifdef SCHEMA_HAS_IfcGradientCurve
+BIND(IfcGradientCurve);


@aothms This is a little problematic. IfcGradientCurve is a IfcCompositeCurve so the item gets processed twice and the result for IfcCompositeCurve is returned. When processed as IfCompositeCurve the IfcCurveSegment are the vertical segments.

I have some code in my driver program that maps a IfcGradientCurve but gets back the implicit_item for the vertical components.

if (item->as<Ifc4x3_add1::IfcGradientCurve>()) { auto mapped_item = ifcopenshell::geometry::taxonomy::cast<ifcopenshell::geometry::taxonomy::implicit_item>(mapping->map(item->as<Ifc4x3_add1::IfcGradientCurve>())); auto loop = ifcopenshell::geometry::taxonomy::cast<ifcopenshell::geometry::taxonomy::loop>(mapped_item->evaluate()); ...

Reversing the order of BIND(IfcGradientCurve) and BIND(IfcCompositeCurve) fixes the problem, but I don't think it is the best solution.

Oh sorry. Yes, that could be addressed like this: https://github.com/IfcOpenShell/IfcOpenShell/pull/3827/files#r1347391406

It's something I introduced in the implementation of caching. I could no longer use the return statement anymore (because writing to cache is at the bottom of the function). So now it "falls through" to the other. That's why my suggestion of the addition of !item && helps circumvent that (well at least in the cases the mapping is successful in assigning something to item).

@aothms This is still a little problematic. If BIND(IfcCompositeCurve) is listed before BIND(IfcGradientCurve) the IfcGradientCurve is treated as an IfcCompositeCurve and the IfcGradientCurve code never gets called.

Is there a function similar to inst->as<IfcSchema::T>() that tests if the most derived type of an instance is equal to a specific type?

For IfcGradientCurve the test for IfcCompositeCurve would be false since the most derived type is IfcGradientCurve.

Listing the BIND statements for the more derived types before the parent types solves the problem, but this seems like a brittle and problematic solution.

Listing the BIND statements for the more derived types before the parent types solves the problem, but this seems like a brittle and problematic solution.

I agree, but in quite a few cases we do depend on matching also inherited types:

IfcOpenShell/src/ifcgeom/mapping/mapping.i

Lines 33 to 37 in 661c8bc

// IfcFacetedBrep included

// IfcAdvancedBrep included

// IfcFacetedBrepWithVoids included

// IfcAdvancedBrepWithVoids included

BIND(IfcManifoldSolidBrep);

IfcOpenShell/src/ifcgeom/mapping/mapping.i

Lines 103 to 105 in 661c8bc

// IfcFaceSurface included

// IfcAdvancedFace included in case of IFC4

BIND(IfcFace);

On the other hand in many other cases it's also problematic because in the inheritance hierarchy of the newly added IFC geometry types, they don't seem to consider liskov-substitution in determining whether entities are supposed to be subtypes, but instead rather pragmatically try to reuse some attribute definitions.

So maybe we need to enhance the BIND() macro with an option to selectively specify whether inheritance is considered to have the best of both worlds?

We could check equality of &inst->declaration() == &Ifc2x3::IfcWall::Class(). This disregards inheritance (and is also more efficient because it doesn't rely on RTTI and virtual inheritance).

aothms · 2023-10-05T13:05:34Z

src/ifcgeom/mapping/bind_convert_impl.i

@@ -5,7 +5,7 @@
 #define BIND(T) \
 	if (inst->as<IfcSchema::T>()) { \


Suggested change

if (inst->as<IfcSchema::T>()) { \

if (!item && inst->as<IfcSchema::T>()) { \

First attempt at taxonomy::implicit_item

b51b451

aothms requested a review from civilx64 September 29, 2023 20:26

RickBrice reviewed Sep 30, 2023

View reviewed changes

civilx64 reviewed Sep 30, 2023

View reviewed changes

aothms added 5 commits October 1, 2023 10:44

Structured binding

5da9eea

Correct inheritance

622a21c

Correct vertical IfcLine handling

dbea6d0

Use boost::math::constants::pi

906f824

Caching of taxonomy::item in map()

c3c6402

RickBrice reviewed Oct 2, 2023

View reviewed changes

RickBrice reviewed Oct 5, 2023

View reviewed changes

aothms commented Oct 5, 2023

View reviewed changes

Krande mentioned this pull request Oct 12, 2023

clang-format: add styles and reformat ifcparse files #3761

Merged

	// IfcFacetedBrep included
	// IfcAdvancedBrep included
	// IfcFacetedBrepWithVoids included
	// IfcAdvancedBrepWithVoids included
	BIND(IfcManifoldSolidBrep);

	// IfcFaceSurface included
	// IfcAdvancedFace included in case of IFC4
	BIND(IfcFace);

		@@ -5,7 +5,7 @@
		#define BIND(T) \
		if (inst->as<IfcSchema::T>()) { \

	if (inst->as<IfcSchema::T>()) { \
	if (!item && inst->as<IfcSchema::T>()) { \

First attempt at taxonomy::implicit_item #3827

Are you sure you want to change the base?

First attempt at taxonomy::implicit_item #3827

Conversation

aothms commented Sep 29, 2023

RickBrice left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

civilx64 Sep 30, 2023 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

RickBrice Oct 1, 2023 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

civilx64 left a comment

Choose a reason for hiding this comment

aothms commented Oct 1, 2023

RickBrice commented Oct 1, 2023

aothms commented Oct 2, 2023

RickBrice commented Oct 2, 2023

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

civilx64 commented Oct 3, 2023

aothms commented Oct 5, 2023

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

civilx64 Sep 30, 2023 •

edited

RickBrice Oct 1, 2023 •

edited