ALS016 Industry best practice for horizontal alignment geometric cont…

…inuity (IVS-20) (#243) * initial wip checkpoint for ALS016 implementation * revise for row-major ordering of transform returned from ifcos * fix test data * remove curvature considerations for rule version 1; remove duplicate code in step implementation * Update features/steps/thens/alignment.py per review comment Co-authored-by: Thomas Krijnen <[email protected]> * Break "pairwise" into a separate step Co-authored-by: Thomas Krijnen <[email protected]> * add clarity on toloerances used to determine continuity * address review comments. New step implementation needs additional work to limit amount of tuple unpacking * fix for failing tests; address additional review comments * - Edit unit test files to only fail a single scenario - Edit geometry calculations to accomodate multiple types of placements (previously they assumed IfcAxis2Placement2D or IfcAxis2Placement3D, which wouldn't accomodate other valid placement types such as IfcAxis2PlacementLinear) --------- Co-authored-by: Thomas Krijnen <[email protected]>
buildingSMART · Sep 10, 2024 · aef5e4a · aef5e4a
1 parent 9a0429c
commit aef5e4a
Show file tree

Hide file tree

Showing 10 changed files with 1,616 additions and 7 deletions.
diff --git a/features/ALS016_Alignment-horizontal-segment-geometric-continuity.feature b/features/ALS016_Alignment-horizontal-segment-geometric-continuity.feature
@@ -0,0 +1,33 @@
+@industry-practice
+@ALS
+@version1
+
+Feature: ALS016 - Alignment horizontal segment geometric continuity
+
+  The rule verifies that there is geometric continuity between segments in an IfcCompositeCurve.
+  The calculated end position and tangent vector of segment `n` is compared to the provided placement of segment `n + 1`.
+  A warning is emitted if the calculated difference is greater than the applicable tolerance.
+  The tolerance for positional continuity is taken from the precision of the applicable geometric context.
+  The tolerance for tangential continuity is taken from the precision of the applicable geometric context and
+  adjusted based on the length of the alignment segment.
+
+Background:
+
+  Given A model with Schema "IFC4.3"
+  Given An IfcAlignment
+  Given Its attribute Representation
+  Given Its attribute Representations
+  Given RepresentationType = 'Curve2D'
+  Given All referenced instances
+  Given Its Entity Type is 'IfcCompositeCurve'
+  Given Its attribute Segments
+  Given Its Entity Type is 'IfcCurveSegment'
+  Given The values grouped pairwise at depth 1
+
+Scenario: Geometric continuity in position
+
+  Then Each segment must have geometric continuity in position
+
+Scenario: Geometric continuity in tangency
+
+  Then Each segment must have geometric continuity in tangency
diff --git a/features/steps/givens/attributes.py b/features/steps/givens/attributes.py
@@ -1,4 +1,5 @@
 import ast
+import itertools
 import operator
 
 import ifcopenshell
@@ -136,7 +137,8 @@ def step_impl(context, file_or_model, field, values):
 
 @gherkin_ifc.step('Its attribute {attribute}')
 def step_impl(context, inst, attribute, tail="single"):
-    yield ValidationOutcome(instance_id=getattr(inst, attribute, None), severity = OutcomeSeverity.PASSED)
+    yield ValidationOutcome(instance_id=getattr(inst, attribute, None), severity=OutcomeSeverity.PASSED)
+
 
 @gherkin_ifc.step("Its {attribute} attribute {condition} with {prefix}")
 def step_impl(context, inst, attribute, condition, prefix):

diff --git a/features/steps/givens/values.py b/features/steps/givens/values.py
@@ -1,7 +1,17 @@
+import itertools
+
 from validation_handling import gherkin_ifc
 from . import ValidationOutcome, OutcomeSeverity
 
+
 @gherkin_ifc.step("Its values")
 @gherkin_ifc.step("Its values excluding {excluding}")
-def step_impl(context, inst, excluding = None):
-    yield ValidationOutcome(instance_id=inst.get_info(recursive=True, include_identifier=False, ignore=excluding), severity=OutcomeSeverity.PASSED)
+def step_impl(context, inst, excluding=None):
+    yield ValidationOutcome(instance_id=inst.get_info(recursive=True, include_identifier=False, ignore=excluding),
+                            severity=OutcomeSeverity.PASSED)
+
+
+@gherkin_ifc.step("The values grouped pairwise at depth 1")
+def step_impl(context, inst):
+    inst = itertools.pairwise(inst)
+    yield ValidationOutcome(instance_id=inst, severity=OutcomeSeverity.PASSED)
diff --git a/features/steps/thens/alignment.py b/features/steps/thens/alignment.py
@@ -1,8 +1,14 @@
+import math
+
 from behave import register_type
 from typing import Dict, List, Union, Optional
+
 import ifcopenshell.entity_instance
+import ifcopenshell.util.unit
 
 from utils import ifc43x_alignment_validation as ifc43
+from utils import geometry
+from utils import ifc
 from validation_handling import gherkin_ifc
 from . import ValidationOutcome, OutcomeSeverity
 
@@ -170,6 +176,7 @@ def ala003_activation_inst(inst, context) -> Union[ifcopenshell.entity_instance
                     if item.id() == inst.id():
                         return candidate.ShapeOfProduct[0]
 
+
 @gherkin_ifc.step(
     'A representation by {ifc_rep_criteria} requires the {existence:absence_or_presence} of {entities} in the business logic')
 def step_impl(context, inst, ifc_rep_criteria, existence, entities):
@@ -279,9 +286,9 @@ def step_impl(context, inst, activation_phrase):
 
             if activation_ent.is_a().upper() == "IFCALIGNMENT":
                 # ensure that all three representation types will be validated
-                  if inst.is_a().upper() in ["IFCSEGMENTEDREFERENCECURVE", "IFCGRADIENTCURVE"]:
+                if inst.is_a().upper() in ["IFCSEGMENTEDREFERENCECURVE", "IFCGRADIENTCURVE"]:
                     inst = inst.BaseCurve
-            
+
             match activation_phrase:
                 case "segment in the applicable IfcAlignment layout":
                     align = ifc43.entities.Alignment().from_entity(activation_ent)
@@ -352,3 +359,61 @@ def step_impl(context, inst, activation_phrase):
                         observed=observed_msg,
                         severity=OutcomeSeverity.ERROR,
                     )
+
+
+@gherkin_ifc.step('Each segment must have geometric continuity in {continuity_type}')
+def step_impl(context, inst, continuity_type):
+    """
+    Assess geometric continuity between alignment segments for ALS016, ALS017, and ALS018
+    """
+    # Ref: https://standards.buildingsmart.org/IFC/RELEASE/IFC4_3/HTML/lexical/IfcTransitionCode.htm
+    position_transition_codes = ["CONTINUOUS", "CONTSAMEGRADIENT", "CONTSAMEGRADIENTSAMECURVATURE"]
+    tangency_transition_codes = ["CONTSAMEGRADIENT", "CONTSAMEGRADIENTSAMECURVATURE"]
+
+    length_unit_scale_factor = ifcopenshell.util.unit.calculate_unit_scale(
+        ifc_file=context.model,
+        unit_type="LENGTHUNIT"
+    )
+    for previous, current in inst:
+        entity_contexts = ifc.recurrently_get_entity_attr(context, current, 'IfcRepresentation', 'ContextOfItems')
+        precision = ifc.get_precision_from_contexts(entity_contexts)
+
+        # calculate number of significant figures to display
+        # use the precision of the geometric context plus one additional digit to accommodate rounding
+        from math import ceil, log10
+        display_sig_figs = abs(int(ceil(log10(precision)))) + 1
+
+        continuity_type = continuity_type.lower()
+
+        if (continuity_type == "position") and (current.Transition in position_transition_codes):
+            expected = precision
+            observed = geometry.alignment_segment_positional_difference(
+                length_unit_scale_factor,
+                previous,
+                current
+            )
+
+        elif (continuity_type == "tangency") and (current.Transition in tangency_transition_codes):
+            expected = math.atan2(precision, current.SegmentLength.wrappedValue)
+            observed = geometry.alignment_segment_angular_difference(
+                length_unit_scale_factor,
+                previous,
+                current
+            )
+        else:
+            return
+
+        if abs(observed) > expected:
+            yield ValidationOutcome(
+                inst=current,
+                expected={
+                    "expected": expected,
+                    "num_digits": display_sig_figs,
+                    "context": f"max deviation in {continuity_type.lower()}",
+                },
+                observed={
+                    "observed": observed,
+                    "num_digits": display_sig_figs,
+                    "context": f"calculated deviation in {continuity_type.lower()}",
+                },
+                severity=OutcomeSeverity.WARNING)
diff --git a/features/steps/utils/geometry.py b/features/steps/utils/geometry.py
@@ -1,11 +1,18 @@
 import operator
 import math
 
+import numpy as np
+
+import ifcopenshell.entity_instance
+import ifcopenshell.geom as ifcos_geom
+import ifcopenshell.ifcopenshell_wrapper as wrapper
+
 from .misc import is_a
 from .ifc import get_precision_from_contexts, recurrently_get_entity_attr
 
 GEOM_TOLERANCE = 1E-12
 
+
 def get_edges(file, inst, sequence_type=frozenset, oriented=False):
     edge_type = tuple if oriented else frozenset
 
@@ -26,12 +33,12 @@ def inner():
                 for ed in bnd.Bound.EdgeList:
                     # @todo take into account edge geometry
                     # edge_geom = ed[2].EdgeGeometry.get_info(recursive=True, include_identifier=False)
-                    
+
                     coords = [
                         ed.EdgeElement.EdgeStart.VertexGeometry.Coordinates,
                         ed.EdgeElement.EdgeEnd.VertexGeometry.Coordinates,
                     ]
-                    
+
                     # @todo verify:
                     # @tfk: afaict, sense only affects the parametric space of the underlying curve,
                     #       not the topology of the begin/end vertices
@@ -70,6 +77,7 @@ def emit(loop):
                         yield from emit(inner)
         else:
             raise NotImplementedError(f"get_edges({inst.is_a()})")
+
     return sequence_type(inner())
 
 
@@ -95,3 +103,78 @@ def is_closed(context, instance):
     precision = get_precision_from_contexts(entity_contexts)
     points_coordinates = get_points(instance)
     return math.dist(points_coordinates[0], points_coordinates[-1]) < precision
+
+
+def evaluate_segment(segment: ifcopenshell.entity_instance, dist_along: float) -> np.ndarray:
+    """
+    Use ifcopenshell to calculate the 4x4 geometric transform at a point on an alignment segment
+    :param segment: The segment containing the point that we would like to
+    :param dist_along: The distance along this segment at the point of interest (point to be calculated)
+    """
+    s = ifcos_geom.settings()
+    pwf = wrapper.map_shape(s, segment.wrapped_data)
+
+    prev_trans_matrix = pwf.evaluate(dist_along)
+
+    return np.array(prev_trans_matrix, dtype=np.float64).T
+
+
+def alignment_segment_positional_difference(
+        length_unit_scale_factor: float, previous_segment: ifcopenshell.entity_instance,
+        segment_to_analyze: ifcopenshell.entity_instance):
+    """
+    Use ifcopenshell to determine the difference in cartesian position between segments of an IfcAlignment.
+    The expected entity type is either `IfcCurveSegment` or `IfcCompositeCurveSegment`.
+
+    :param length_unit_scale_factor: Scale factor between the project units and metric units used internally by
+    ifcopenshell
+    :param previous_segment: The segment that precede the segment being analyzed.  The end point of this segment
+    will be determined via ifcopenshell geometry calculations.
+    :param segment_to_analyze: The segment under analysis.  The calculated end point of the previous segment will be
+    compared to the calculated start point of this segment.
+    """
+
+    u = abs(previous_segment.SegmentLength.wrappedValue) * length_unit_scale_factor
+    prev_end_transform = evaluate_segment(segment=previous_segment, dist_along=u)
+    current_start_transform = evaluate_segment(segment=segment_to_analyze, dist_along=0.0)
+
+    e0 = prev_end_transform[3][0] / length_unit_scale_factor
+    e1 = prev_end_transform[3][1] / length_unit_scale_factor
+    preceding_end = (e0, e1)
+
+    s0 = current_start_transform[3][0] / length_unit_scale_factor
+    s1 = current_start_transform[3][1] / length_unit_scale_factor
+    current_start = (s0, s1)
+
+    return math.dist(preceding_end, current_start)
+
+
+def alignment_segment_angular_difference(
+        length_unit_scale_factor: float, previous_segment: ifcopenshell.entity_instance,
+        segment_to_analyze: ifcopenshell.entity_instance):
+    """
+    Use ifcopenshell to determine the difference in tangent direction angle between segments of an IfcAlignment.
+    The expected entity type is either `IfcCurveSegment` or `IfcCompositeCurveSegment`.
+
+    :param length_unit_scale_factor: Scale factor between the project units and metric units used internally by
+    ifcopenshell
+    :param previous_segment: The segment that precede the segment being analyzed.  The ending direction of this segment
+    will be determined via ifcopenshell geometry calculations.
+    :param segment_to_analyze: The segment under analysis.  The calculated ending direction of the previous segment
+    will be compared to the calculated starting direction of this segment.
+    """
+    u = abs(float(previous_segment.SegmentLength.wrappedValue)) * length_unit_scale_factor
+    prev_end_transform = evaluate_segment(segment=previous_segment, dist_along=u)
+    current_start_transform = evaluate_segment(segment=segment_to_analyze, dist_along=0.0)
+
+    prev_i = prev_end_transform[0][0]
+    prev_j = prev_end_transform[0][1]
+    preceding_end_direction = math.atan2(prev_j, prev_i)
+
+    curr_i = current_start_transform[0][0]
+    curr_j = current_start_transform[0][1]
+    current_start_direction = math.atan2(curr_j, curr_i)
+
+    delta = abs(current_start_direction - preceding_end_direction)
+
+    return delta