"""Compare the occupied 3D geometry of two STEP files with OpenCASCADE.""" from __future__ import annotations import itertools import math from pathlib import Path from typing import Any import numpy as np from OCP.Bnd import Bnd_Box from OCP.BRep import BRep_Builder, BRep_Tool from OCP.BRepAlgoAPI import BRepAlgoAPI_Common from OCP.BRepBndLib import BRepBndLib from OCP.BRepBuilderAPI import BRepBuilderAPI_MakeVertex, BRepBuilderAPI_Transform from OCP.BRepCheck import BRepCheck_Analyzer from OCP.BRepExtrema import BRepExtrema_DistShapeShape from OCP.BRepGProp import BRepGProp from OCP.BRepMesh import BRepMesh_IncrementalMesh from OCP.gp import gp_Pnt, gp_Trsf from OCP.GProp import GProp_GProps from OCP.IFSelect import IFSelect_RetDone from OCP.STEPControl import STEPControl_Reader from OCP.TopAbs import ( TopAbs_EDGE, TopAbs_FACE, TopAbs_SHELL, TopAbs_SOLID, TopAbs_VERTEX, ) from OCP.TopExp import TopExp_Explorer from OCP.TopLoc import TopLoc_Location from OCP.TopoDS import TopoDS, TopoDS_Compound, TopoDS_Shape LENGTH_TOLERANCE_MM = 1.0e-2 VOLUME_ABSOLUTE_TOLERANCE_MM3 = 1.0 VOLUME_RELATIVE_TOLERANCE = 1.0e-5 RESIDUAL_VOLUME_ABSOLUTE_TOLERANCE_MM3 = 1.0e-1 AREA_ABSOLUTE_TOLERANCE_MM2 = 1.0 AREA_RELATIVE_TOLERANCE = 1.0e-5 INERTIA_RELATIVE_TOLERANCE = 1.0e-5 def read_step(path: Path) -> TopoDS_Shape: reader = STEPControl_Reader() if reader.ReadFile(str(path)) != IFSelect_RetDone: raise RuntimeError(f"cannot read STEP file: {path}") if reader.TransferRoots() < 1: raise RuntimeError(f"STEP file has no transferable roots: {path}") shape = reader.OneShape() if shape.IsNull(): raise RuntimeError(f"STEP import produced a null shape: {path}") return shape def topology_count(shape: TopoDS_Shape, shape_type: Any) -> int: explorer = TopExp_Explorer(shape, shape_type) count = 0 while explorer.More(): count += 1 explorer.Next() return count def downcast(shape: TopoDS_Shape, caster_name: str) -> TopoDS_Shape: caster = getattr(TopoDS, caster_name, None) if caster is None: caster = getattr(TopoDS, f"{caster_name}_s") return caster(shape) def only_solid(shape: TopoDS_Shape) -> TopoDS_Shape | None: explorer = TopExp_Explorer(shape, TopAbs_SOLID) if not explorer.More(): return None solid = downcast(explorer.Current(), "Solid") explorer.Next() return None if explorer.More() else solid def maximum_topological_tolerance(shape: TopoDS_Shape) -> float: maximum = 0.0 for shape_type, caster_name in ( (TopAbs_VERTEX, "Vertex"), (TopAbs_EDGE, "Edge"), (TopAbs_FACE, "Face"), ): # OCP 7 uses Face_s/Edge_s/Vertex_s; OCP 8 removed the suffix. explorer = TopExp_Explorer(shape, shape_type) while explorer.More(): maximum = max( maximum, BRep_Tool.Tolerance_s(downcast(explorer.Current(), caster_name)), ) explorer.Next() return maximum def metrics(shape: TopoDS_Shape) -> dict[str, Any]: volume_properties = GProp_GProps() volume_error = BRepGProp.VolumeProperties_s(shape, volume_properties, 1.0e-10, True, False) centre = volume_properties.CentreOfMass() # Recompute inertia after translating the body close to the origin. This # avoids catastrophic cancellation for otherwise-identical STEP files # whose author chose coordinates far from (0, 0, 0). centering = gp_Trsf() centering.SetValues( 1.0, 0.0, 0.0, -centre.X(), 0.0, 1.0, 0.0, -centre.Y(), 0.0, 0.0, 1.0, -centre.Z(), ) centred_shape = BRepBuilderAPI_Transform(shape, centering, True, False).Shape() centred_volume_properties = GProp_GProps() BRepGProp.VolumeProperties_s(centred_shape, centred_volume_properties, 1.0e-10, True, False) surface_properties = GProp_GProps() surface_error = BRepGProp.SurfaceProperties_s(shape, surface_properties, 1.0e-10, False) bounds = Bnd_Box() bounds.SetGap(0.0) BRepBndLib.AddOptimal_s(shape, bounds, True, False) inertia = centred_volume_properties.MatrixOfInertia() inertia_matrix = np.array( [[inertia.Value(row, column) for column in range(1, 4)] for row in range(1, 4)], dtype=float, ) principal_moments, principal_axes = np.linalg.eigh(inertia_matrix) if hasattr(bounds, "GetXMin"): bbox = [ bounds.GetXMin(), bounds.GetYMin(), bounds.GetZMin(), bounds.GetXMax(), bounds.GetYMax(), bounds.GetZMax(), ] else: bbox = list(bounds.Get()) solid_count = topology_count(shape, TopAbs_SOLID) shell_count = topology_count(shape, TopAbs_SHELL) face_count = topology_count(shape, TopAbs_FACE) edge_count = topology_count(shape, TopAbs_EDGE) vertex_count = topology_count(shape, TopAbs_VERTEX) solid = only_solid(shape) solid_face_count = topology_count(solid, TopAbs_FACE) if solid is not None else 0 solid_edge_count = topology_count(solid, TopAbs_EDGE) if solid is not None else 0 solid_vertex_count = topology_count(solid, TopAbs_VERTEX) if solid is not None else 0 return { "valid_brep": bool(BRepCheck_Analyzer(shape, True).IsValid()), "solid_count": solid_count, "shell_count": shell_count, "face_count": face_count, "extra_face_count_outside_only_solid": max(0, face_count - solid_face_count), "extra_edge_count_outside_only_solid": max(0, edge_count - solid_edge_count), "extra_vertex_count_outside_only_solid": max(0, vertex_count - solid_vertex_count), "maximum_topological_tolerance_mm": maximum_topological_tolerance(shape), "volume_mm3": volume_properties.Mass(), "surface_area_mm2": surface_properties.Mass(), "centre_of_mass_mm": [centre.X(), centre.Y(), centre.Z()], "inertia_matrix_at_com_mm5": inertia_matrix.tolist(), "principal_moments_mm5": principal_moments.tolist(), "principal_axes": principal_axes.tolist(), "bounding_box_mm": bbox, "estimated_relative_integration_error": { "volume": volume_error, "surface_area": surface_error, }, } def intersection( left: TopoDS_Shape, right: TopoDS_Shape, fuzzy_tolerance: float, label: str, ) -> TopoDS_Shape: operation = BRepAlgoAPI_Common(left, right) operation.SetNonDestructive(True) operation.SetFuzzyValue(fuzzy_tolerance) operation.Build() if not operation.IsDone(): raise RuntimeError(f"OpenCASCADE Boolean intersection failed: {label}") return operation.Shape() def volume(shape: TopoDS_Shape) -> float: if shape.IsNull() or topology_count(shape, TopAbs_SOLID) == 0: return 0.0 properties = GProp_GProps() BRepGProp.VolumeProperties_s(shape, properties, 1.0e-10, True, False) return abs(properties.Mass()) def principal_moments_match( metrics_a: dict[str, Any], metrics_b: dict[str, Any], relative_tolerance: float, ) -> bool: moments_a = np.asarray(metrics_a["principal_moments_mm5"]) moments_b = np.asarray(metrics_b["principal_moments_mm5"]) scale = max( float(np.max(np.abs(moments_a))), float(np.max(np.abs(moments_b))), 1.0, ) return bool(np.max(np.abs(moments_a - moments_b)) <= relative_tolerance * scale) def inertia_is_degenerate(model_metrics: dict[str, Any], relative_tolerance: float) -> bool: moments = np.sort(np.asarray(model_metrics["principal_moments_mm5"])) scale = max(float(np.max(np.abs(moments))), 1.0) return bool(np.any(np.diff(moments) <= relative_tolerance * scale)) def rigid_alignment_candidates( metrics_a: dict[str, Any], metrics_b: dict[str, Any], inertia_relative_tolerance: float, ) -> list[dict[str, Any]]: """Generate all proper rotations consistent with the principal inertia axes.""" centre_a = np.asarray(metrics_a["centre_of_mass_mm"], dtype=float) centre_b = np.asarray(metrics_b["centre_of_mass_mm"], dtype=float) moments_a = np.asarray(metrics_a["principal_moments_mm5"], dtype=float) moments_b = np.asarray(metrics_b["principal_moments_mm5"], dtype=float) axes_a = np.asarray(metrics_a["principal_axes"], dtype=float) axes_b = np.asarray(metrics_b["principal_axes"], dtype=float) moment_scale = max( float(np.max(np.abs(moments_a))), float(np.max(np.abs(moments_b))), 1.0, ) moment_tolerance = inertia_relative_tolerance * moment_scale candidates: list[dict[str, Any]] = [] for permutation in itertools.permutations(range(3)): mapped_moments = moments_b[list(permutation)] if np.max(np.abs(moments_a - mapped_moments)) > moment_tolerance: continue for signs in itertools.product((-1.0, 1.0), repeat=3): signed_permutation = np.zeros((3, 3), dtype=float) for row, column in enumerate(permutation): signed_permutation[row, column] = signs[row] rotation = axes_a @ signed_permutation @ axes_b.T # Only SO(3) is allowed: reflections and scaling are forbidden. if np.linalg.det(rotation) < 0.999999999: continue translation = centre_a - rotation @ centre_b if any( np.max(np.abs(rotation - existing["rotation_array"])) < 1.0e-10 for existing in candidates ): continue candidates.append( { "rotation_array": rotation, "translation_array": translation, "source": "volume principal inertia axes", } ) return candidates def apply_rigid_transform( shape: TopoDS_Shape, rotation: np.ndarray, translation: np.ndarray ) -> TopoDS_Shape: transform = gp_Trsf() transform.SetValues( float(rotation[0, 0]), float(rotation[0, 1]), float(rotation[0, 2]), float(translation[0]), float(rotation[1, 0]), float(rotation[1, 1]), float(rotation[1, 2]), float(translation[1]), float(rotation[2, 0]), float(rotation[2, 1]), float(rotation[2, 2]), float(translation[2]), ) if transform.IsNegative() or abs(transform.ScaleFactor() - 1.0) > 1.0e-10: raise RuntimeError("candidate transform is not a proper rigid transform") result = BRepBuilderAPI_Transform(shape, transform, True, False).Shape() if result.IsNull(): raise RuntimeError("OpenCASCADE rigid transformation returned a null shape") return result def evenly_spaced_indices(count: int, limit: int) -> list[int]: """Return deterministic 1-based indices spanning an OpenCASCADE array.""" if count <= limit: return list(range(1, count + 1)) return sorted({1 + round(position * (count - 1) / (limit - 1)) for position in range(limit)}) def boundary_shape(shape: TopoDS_Shape) -> TopoDS_Compound: """Build a compound containing only the faces (the solid boundary).""" boundary = TopoDS_Compound() builder = BRep_Builder() builder.MakeCompound(boundary) explorer = TopExp_Explorer(shape, TopAbs_FACE) while explorer.More(): builder.Add(boundary, explorer.Current()) explorer.Next() return boundary def boundary_sample_points( shape: TopoDS_Shape, mesh_deflection: float, max_nodes_per_face: int = 2048, triangle_samples_per_face: int = 24, ) -> list[gp_Pnt]: """Sample each B-rep face with mesh nodes and triangle interior points.""" mesher = BRepMesh_IncrementalMesh(shape, mesh_deflection, False, 0.05, False) mesher.Perform() if not mesher.IsDone(): raise RuntimeError("OpenCASCADE boundary meshing failed") points: list[gp_Pnt] = [] seen: set[tuple[float, float, float]] = set() explorer = TopExp_Explorer(shape, TopAbs_FACE) while explorer.More(): face = downcast(explorer.Current(), "Face") location = TopLoc_Location() triangulation = BRep_Tool.Triangulation_s(face, location) if triangulation is None or triangulation.NbNodes() < 1: raise RuntimeError("OpenCASCADE produced a face without triangulation") location_transform = location.Transformation() def add_point(point: gp_Pnt, transform=location_transform) -> None: transformed = point.Transformed(transform) key = ( round(transformed.X(), 10), round(transformed.Y(), 10), round(transformed.Z(), 10), ) if key not in seen: seen.add(key) points.append(transformed) for node_index in evenly_spaced_indices(triangulation.NbNodes(), max_nodes_per_face): add_point(triangulation.Node(node_index)) for triangle_index in evenly_spaced_indices( triangulation.NbTriangles(), triangle_samples_per_face ): node_indices = triangulation.Triangle(triangle_index).Get() nodes = [triangulation.Node(index) for index in node_indices] add_point( gp_Pnt( sum(point.X() for point in nodes) / 3.0, sum(point.Y() for point in nodes) / 3.0, sum(point.Z() for point in nodes) / 3.0, ) ) explorer.Next() if not points: raise RuntimeError("OpenCASCADE produced no boundary sample points") return points def directed_boundary_sample_distance( source: TopoDS_Shape, target: TopoDS_Shape, length_tolerance: float, ) -> dict[str, Any]: """Measure sampled source-boundary points against the exact target boundary.""" mesh_deflection = max(1.0e-7, length_tolerance * 0.25) points = boundary_sample_points(source, mesh_deflection) target_boundary = boundary_shape(target) extrema = BRepExtrema_DistShapeShape() extrema.LoadS2(target_boundary) extrema.SetDeflection(max(1.0e-10, length_tolerance * 1.0e-4)) observed_maximum = 0.0 evaluated = 0 for point in points: vertex = BRepBuilderAPI_MakeVertex(point).Vertex() extrema.LoadS1(vertex) extrema.Perform() if not extrema.IsDone(): raise RuntimeError("OpenCASCADE point-to-boundary distance failed") distance = float(extrema.Value()) if not math.isfinite(distance): raise RuntimeError("OpenCASCADE returned a non-finite boundary distance") observed_maximum = max(observed_maximum, distance) evaluated += 1 if observed_maximum > length_tolerance: break return { "within_tolerance": observed_maximum <= length_tolerance, "observed_maximum_mm": observed_maximum, "evaluated_sample_count": evaluated, "available_sample_count": len(points), "mesh_deflection_mm": mesh_deflection, "early_exit": evaluated < len(points), } def serializable_transform(candidate: dict[str, Any]) -> dict[str, Any]: rotation = candidate["rotation_array"] translation = candidate["translation_array"] matrix = np.eye(4) matrix[:3, :3] = rotation matrix[:3, 3] = translation return { "source": candidate["source"], "rotation_3x3": rotation.tolist(), "translation_mm": translation.tolist(), "homogeneous_4x4": matrix.tolist(), "rotation_determinant": float(np.linalg.det(rotation)), "orthogonality_error": float(np.max(np.abs(rotation.T @ rotation - np.eye(3)))), } def compare( path_a: Path, path_b: Path, ) -> dict[str, Any]: shape_a = read_step(path_a) shape_b = read_step(path_b) metrics_a = metrics(shape_a) metrics_b = metrics(shape_b) effective_volume_tolerance = max( VOLUME_ABSOLUTE_TOLERANCE_MM3, VOLUME_RELATIVE_TOLERANCE * max(abs(metrics_a["volume_mm3"]), abs(metrics_b["volume_mm3"]), 1.0), ) effective_area_tolerance = max( AREA_ABSOLUTE_TOLERANCE_MM2, AREA_RELATIVE_TOLERANCE * max( abs(metrics_a["surface_area_mm2"]), abs(metrics_b["surface_area_mm2"]), 1.0, ), ) invariant_checks = { "a_is_valid_solid": metrics_a["valid_brep"] and metrics_a["solid_count"] == 1 and metrics_a["extra_face_count_outside_only_solid"] == 0 and metrics_a["extra_edge_count_outside_only_solid"] == 0 and metrics_a["extra_vertex_count_outside_only_solid"] == 0 and math.isfinite(metrics_a["volume_mm3"]) and metrics_a["volume_mm3"] > 0.0, "b_is_valid_solid": metrics_b["valid_brep"] and metrics_b["solid_count"] == 1 and metrics_b["extra_face_count_outside_only_solid"] == 0 and metrics_b["extra_edge_count_outside_only_solid"] == 0 and metrics_b["extra_vertex_count_outside_only_solid"] == 0 and math.isfinite(metrics_b["volume_mm3"]) and metrics_b["volume_mm3"] > 0.0, "volumes_match": abs(metrics_a["volume_mm3"] - metrics_b["volume_mm3"]) <= effective_volume_tolerance, "surface_areas_match": abs(metrics_a["surface_area_mm2"] - metrics_b["surface_area_mm2"]) <= effective_area_tolerance, "principal_moments_match": principal_moments_match( metrics_a, metrics_b, INERTIA_RELATIVE_TOLERANCE ), "shell_counts_match": metrics_a["shell_count"] == metrics_b["shell_count"], } degenerate_inertia = inertia_is_degenerate( metrics_a, INERTIA_RELATIVE_TOLERANCE ) or inertia_is_degenerate(metrics_b, INERTIA_RELATIVE_TOLERANCE) candidates = ( rigid_alignment_candidates(metrics_a, metrics_b, INERTIA_RELATIVE_TOLERANCE) if invariant_checks["principal_moments_match"] else [] ) input_model_tolerance = max( metrics_a["maximum_topological_tolerance_mm"], metrics_b["maximum_topological_tolerance_mm"], ) # Cross-CAD STEP files often carry edge tolerances larger than the kernel's # default confusion tolerance. A Boolean tolerance below the input model # tolerance can misclassify almost-coincident solids as completely disjoint. boolean_fuzzy_tolerance = max( 1.0e-10, input_model_tolerance, LENGTH_TOLERANCE_MM * 0.01, ) input_tolerance_compatible = input_model_tolerance <= LENGTH_TOLERANCE_MM posed_candidates: list[dict[str, Any]] = [] for index, candidate in enumerate(candidates): aligned_shape = apply_rigid_transform( shape_b, candidate["rotation_array"], candidate["translation_array"], ) aligned_metrics = metrics(aligned_shape) posed_candidates.append( { "index": index, "candidate": candidate, "shape": aligned_shape, "metrics": aligned_metrics, "bbox_error_mm": max( abs(a - b) for a, b in zip( metrics_a["bounding_box_mm"], aligned_metrics["bounding_box_mm"] ) ), } ) posed_candidates.sort(key=lambda item: item["bbox_error_mm"]) bbox_compatible = [ item for item in posed_candidates if item["bbox_error_mm"] <= LENGTH_TOLERANCE_MM ] # If no pose passes the cheap bounding-box test, evaluate the closest one # anyway so the report still includes a meaningful Boolean residual. candidates_to_evaluate = bbox_compatible or posed_candidates[:1] attempts: list[dict[str, Any]] = [] best: dict[str, Any] | None = None boolean_failures = 0 for item in candidates_to_evaluate: try: common_shape = intersection( shape_a, item["shape"], boolean_fuzzy_tolerance, f"A intersect aligned B, candidate {item['index']}", ) common_volume = volume(common_shape) common_has_solid = ( not common_shape.IsNull() and topology_count(common_shape, TopAbs_SOLID) > 0 ) boolean_retry_used = False # OCCT can occasionally return a completed but empty Common for # nearly coincident cross-CAD B-reps. Retry once with a freshly # transformed operand before accepting that result as disjoint. if ( not common_has_solid and all(invariant_checks.values()) and input_tolerance_compatible and item["bbox_error_mm"] <= LENGTH_TOLERANCE_MM ): boolean_retry_used = True retry_shape = apply_rigid_transform( shape_b, item["candidate"]["rotation_array"], item["candidate"]["translation_array"], ) retry_common_shape = intersection( shape_a, retry_shape, boolean_fuzzy_tolerance, f"A intersect aligned B retry, candidate {item['index']}", ) retry_common_volume = volume(retry_common_shape) if retry_common_volume > common_volume: item["shape"] = retry_shape common_shape = retry_common_shape common_volume = retry_common_volume common_has_solid = ( not common_shape.IsNull() and topology_count(common_shape, TopAbs_SOLID) > 0 ) common_excess = max( 0.0, common_volume - min(metrics_a["volume_mm3"], item["metrics"]["volume_mm3"]), ) if common_excess > RESIDUAL_VOLUME_ABSOLUTE_TOLERANCE_MM3: raise RuntimeError("intersection volume exceeds an input volume beyond tolerance") # For regular solids these identities are exact and are more stable # across exporters than directly cutting nearly coincident faces. a_minus_b = max(0.0, metrics_a["volume_mm3"] - common_volume) b_minus_a = max(0.0, item["metrics"]["volume_mm3"] - common_volume) symmetric_difference = a_minus_b + b_minus_a preliminary_match = ( all(invariant_checks.values()) and input_tolerance_compatible and item["bbox_error_mm"] <= LENGTH_TOLERANCE_MM and a_minus_b <= RESIDUAL_VOLUME_ABSOLUTE_TOLERANCE_MM3 and b_minus_a <= RESIDUAL_VOLUME_ABSOLUTE_TOLERANCE_MM3 ) boundary_a_to_b = None boundary_b_to_a = None if preliminary_match: boundary_a_to_b = directed_boundary_sample_distance( shape_a, item["shape"], LENGTH_TOLERANCE_MM ) boundary_b_to_a = directed_boundary_sample_distance( item["shape"], shape_a, LENGTH_TOLERANCE_MM ) boundary_matches = bool( boundary_a_to_b and boundary_a_to_b["within_tolerance"] and boundary_b_to_a and boundary_b_to_a["within_tolerance"] ) attempt = { "candidate_index": item["index"], "bounding_box_max_error_mm": item["bbox_error_mm"], "a_minus_b_volume_mm3": a_minus_b, "b_minus_a_volume_mm3": b_minus_a, "intersection_volume_mm3": common_volume, "intersection_excess_mm3": common_excess, "symmetric_difference_volume_mm3": symmetric_difference, "sampled_boundary_a_to_b": boundary_a_to_b, "sampled_boundary_b_to_a": boundary_b_to_a, "boolean_succeeded": True, "boolean_retry_used": boolean_retry_used, "intersection_has_solid": common_has_solid, } attempts.append(attempt) evaluated = { **item, "a_minus_b": a_minus_b, "b_minus_a": b_minus_a, "symmetric_difference": symmetric_difference, "boundary_a_to_b": boundary_a_to_b, "boundary_b_to_a": boundary_b_to_a, } if preliminary_match and boundary_matches: best = evaluated break if best is None or symmetric_difference < best["symmetric_difference"]: best = evaluated except RuntimeError as error: boolean_failures += 1 attempts.append( { "candidate_index": item["index"], "bounding_box_max_error_mm": item["bbox_error_mm"], "boolean_succeeded": False, "error": str(error), } ) if best is None: aligned_metrics = None chosen_transform = None a_minus_b = None b_minus_a = None symmetric_difference = None pose_checks = { "bounding_boxes_match_after_alignment": False, "a_minus_b_is_empty_within_tolerance": False, "b_minus_a_is_empty_within_tolerance": False, "sampled_a_boundary_matches_b": False, "sampled_b_boundary_matches_a": False, } else: aligned_metrics = best["metrics"] chosen_transform = serializable_transform(best["candidate"]) a_minus_b = best["a_minus_b"] b_minus_a = best["b_minus_a"] symmetric_difference = best["symmetric_difference"] pose_checks = { "bounding_boxes_match_after_alignment": best["bbox_error_mm"] <= LENGTH_TOLERANCE_MM, "a_minus_b_is_empty_within_tolerance": a_minus_b <= RESIDUAL_VOLUME_ABSOLUTE_TOLERANCE_MM3, "b_minus_a_is_empty_within_tolerance": b_minus_a <= RESIDUAL_VOLUME_ABSOLUTE_TOLERANCE_MM3, "sampled_a_boundary_matches_b": bool( best["boundary_a_to_b"] and best["boundary_a_to_b"]["within_tolerance"] ), "sampled_b_boundary_matches_a": bool( best["boundary_b_to_a"] and best["boundary_b_to_a"]["within_tolerance"] ), } checks = {**invariant_checks, **pose_checks} checks["input_model_tolerances_within_length_tolerance"] = input_tolerance_compatible passed = all(checks.values()) invariant_mismatch = not all(invariant_checks.values()) if passed: status = "SAME" same_geometry: bool | None = True interpretation = "congruent by a proper rigid transform, within tolerance" elif not input_tolerance_compatible: status = "INCONCLUSIVE" same_geometry = None interpretation = ( "an input STEP carries a modeling tolerance larger than the requested length tolerance" ) elif invariant_mismatch: status = "DIFFERENT" same_geometry = False interpretation = "rigid-transform invariants differ" elif boolean_failures == len(candidates_to_evaluate) and candidates_to_evaluate: status = "INCONCLUSIVE" same_geometry = None interpretation = "OpenCASCADE Boolean verification failed" elif degenerate_inertia: status = "INCONCLUSIVE" same_geometry = None interpretation = ( "inertia has repeated principal moments; a geometry landmark or global " "registration fallback is required to exhaust all rotations" ) else: status = "DIFFERENT" same_geometry = False interpretation = "no proper rigid transform produced matching occupied solids" return { "status": status, "same_geometry": same_geometry, "interpretation": interpretation, "alignment_mode": "automatic translation and proper rotation (no reflection or scaling)", "file_a": str(path_a), "file_b": str(path_b), "tolerances": { "length_mm": LENGTH_TOLERANCE_MM, "volume_absolute_mm3": VOLUME_ABSOLUTE_TOLERANCE_MM3, "volume_relative": VOLUME_RELATIVE_TOLERANCE, "effective_volume_mm3": effective_volume_tolerance, "residual_volume_absolute_mm3": RESIDUAL_VOLUME_ABSOLUTE_TOLERANCE_MM3, "area_absolute_mm2": AREA_ABSOLUTE_TOLERANCE_MM2, "area_relative": AREA_RELATIVE_TOLERANCE, "effective_area_mm2": effective_area_tolerance, "boolean_fuzzy_mm": boolean_fuzzy_tolerance, "maximum_input_topological_tolerance_mm": input_model_tolerance, "boundary_mesh_deflection_mm": max(1.0e-7, LENGTH_TOLERANCE_MM * 0.25), "inertia_relative": INERTIA_RELATIVE_TOLERANCE, }, "inertia_degenerate": degenerate_inertia, "alignment": { "candidate_count": len(candidates), "bbox_compatible_candidate_count": len(bbox_compatible), "boolean_evaluated_candidate_count": len(attempts), "chosen_transform_b_to_a": chosen_transform, "attempts": attempts, }, "checks": checks, "differences": { "absolute_volume_mm3": abs(metrics_a["volume_mm3"] - metrics_b["volume_mm3"]), "absolute_surface_area_mm2": abs( metrics_a["surface_area_mm2"] - metrics_b["surface_area_mm2"] ), "a_minus_b_volume_mm3": a_minus_b, "b_minus_a_volume_mm3": b_minus_a, "symmetric_difference_volume_mm3": symmetric_difference, }, "a": metrics_a, "b": metrics_b, "b_after_alignment": aligned_metrics, "note": ( "Face and edge partitioning may differ while same_geometry is true; " "STEP text, entity order, names, colours, and metadata are not compared. " "The linear boundary check samples every face using a fine OpenCASCADE " "triangulation and measures those points against the opposite exact B-rep. " "Automatic alignment never uses reflection or scaling. A failed search on " "a repeated-inertia model is reported as INCONCLUSIVE rather than DIFFERENT." ), }