From eec1aad625f4fdfb0f2945fd9917aafbedcceffc Mon Sep 17 00:00:00 2001 From: Michael Jackson Date: Mon, 13 Jul 2026 13:10:40 -0400 Subject: [PATCH 1/2] BUG: Fix ComputeTriangleGeomCentroids periodic centroid (issue #1665) The periodic path added a constant |max-min|/2 offset to the naive centroid, which is only correct when a wrapped feature's mass is symmetric about the seam and can place the centroid outside the domain for asymmetric features. * Reimplement the AdjustCentroidsForPeriodicFaces BoundingBox overload as a largest-empty-gap minimum-image mean computed from the feature's vertex coordinates, with in-domain wrap-back and an arithmetic-mean fallback for domain-filling features. * Fix a latent std::bad_function_call: the algorithm invoked m_MessageHandler.m_Callback() directly, bypassing the empty-callback guard; use the call operator instead. * Add cancel checks to the triangle and feature loops; avoid a per-feature std::set copy. * Add a Class 1 (Analytical) + Class 4 (Invariant) periodic oracle test covering non-wrapping, symmetric-wrap, asymmetric-wrap, and domain-filling features; passes in-core and out-of-core. * Add V&V deliverables (report, deviations, provenance) and document the Is Periodic behavior. Signed-off-by: Michael Jackson --- .../ComputeTriangleGeomCentroidsFilter.md | 14 ++ .../ComputeTriangleGeomCentroids.cpp | 17 +- .../test/ComputeTriangleGeomCentroidsTest.cpp | 145 ++++++++++++++++++ .../vv/ComputeTriangleGeomCentroidsFilter.md | 84 ++++++++++ .../ComputeTriangleGeomCentroidsFilter.md | 49 ++++++ .../vv/provenance/12_IN625_GBCD.md | 40 +++++ src/simplnx/Utilities/GeometryHelpers.cpp | 131 ++++++++++++++-- src/simplnx/Utilities/GeometryHelpers.hpp | 24 +-- 8 files changed, 479 insertions(+), 25 deletions(-) create mode 100644 src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md create mode 100644 src/Plugins/SimplnxCore/vv/deviations/ComputeTriangleGeomCentroidsFilter.md create mode 100644 src/Plugins/SimplnxCore/vv/provenance/12_IN625_GBCD.md diff --git a/src/Plugins/SimplnxCore/docs/ComputeTriangleGeomCentroidsFilter.md b/src/Plugins/SimplnxCore/docs/ComputeTriangleGeomCentroidsFilter.md index 8928f755cb..2094f23449 100644 --- a/src/Plugins/SimplnxCore/docs/ComputeTriangleGeomCentroidsFilter.md +++ b/src/Plugins/SimplnxCore/docs/ComputeTriangleGeomCentroidsFilter.md @@ -14,6 +14,20 @@ using the following algorithm: node to be entered once for a given owner*) 3. For each **Feature**, find the average (x,y,z) coordinate from the set of nodes that bound it +## Is Periodic + +When **Is Periodic** is enabled, features whose nodes wrap across opposing boundaries of the mesh are handled with a +minimum-image (periodic) centroid rather than the plain average. For any axis on which a feature touches both opposing +faces of the geometry's bounding box, the centroid component on that axis is computed by unwrapping the feature's nodes +across the boundary (using the largest empty gap in their distribution) before averaging, and mapping the result back +into the domain. This depends on where the feature's mass actually lies, so it is correct for asymmetrically-wrapped +features and never places the centroid outside the domain. Axes on which the feature does not wrap keep the ordinary +average, and a feature that fills the whole domain along an axis falls back to the ordinary average on that axis. When +**Is Periodic** is disabled (the default), the plain average from step 3 is reported. + +> **Note:** Enabling **Is Periodic** assumes the periodic domain equals the geometry's bounding box. An informational +> message is emitted for each feature that is adjusted, since a manual review may be warranted. + % Auto generated parameter table will be inserted here ## Example Pipelines diff --git a/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeTriangleGeomCentroids.cpp b/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeTriangleGeomCentroids.cpp index ec08855a0d..3de6c2857a 100644 --- a/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeTriangleGeomCentroids.cpp +++ b/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeTriangleGeomCentroids.cpp @@ -58,6 +58,10 @@ Result<> ComputeTriangleGeomCentroids::operator()() for(MeshIndexType i = 0; i < numTriangles; i++) { + if(m_ShouldCancel) + { + return {}; + } const int32 faceLabel0 = faceLabels[2 * i + 0]; const int32 faceLabel1 = faceLabels[2 * i + 1]; if(faceLabel0 > 0) @@ -78,7 +82,11 @@ Result<> ComputeTriangleGeomCentroids::operator()() for(MeshIndexType i = 0; i < numFeatures; i++) { - std::set vertexSet = vertexSets[i]; + if(m_ShouldCancel) + { + return {}; + } + const std::set& vertexSet = vertexSets[i]; auto periodicFaces = GeometryHelpers::Topology::FindElementPeriodicFaces(boundingBox, vertexCoords, vertexSet); for(const auto& vert : vertexSets[i]) @@ -95,10 +103,11 @@ Result<> ComputeTriangleGeomCentroids::operator()() if(m_InputValues->IsPeriodic) { - if(GeometryHelpers::Topology::AdjustCentroidsForPeriodicFaces(boundingBox, periodicFaces, centroids, i)) + if(GeometryHelpers::Topology::AdjustCentroidsForPeriodicFaces(boundingBox, periodicFaces, vertexCoords, vertexSets[i], centroids, i)) { - IFilter::Message warningMsg{IFilter::Message::Type::Info, fmt::format("Feature ID {} may be periodic. Manual review may be necessary.", i)}; - m_MessageHandler.m_Callback(warningMsg); + // Use the MessageHandler's call operator (which guards against an empty callback) rather than + // invoking m_Callback directly, which throws std::bad_function_call when no handler is installed. + m_MessageHandler(IFilter::Message{IFilter::Message::Type::Info, fmt::format("Feature ID {} may be periodic. Manual review may be necessary.", i)}); } } } diff --git a/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp b/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp index 485945ca1f..0258a4643e 100644 --- a/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp +++ b/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp @@ -1,6 +1,9 @@ #include "SimplnxCore/Filters/ComputeTriangleGeomCentroidsFilter.hpp" #include "SimplnxCore/SimplnxCore_test_dirs.hpp" +#include "simplnx/DataStructure/AttributeMatrix.hpp" +#include "simplnx/DataStructure/DataStore.hpp" +#include "simplnx/DataStructure/Geometry/TriangleGeom.hpp" #include "simplnx/Parameters/ArrayCreationParameter.hpp" #include "simplnx/Parameters/DataObjectNameParameter.hpp" #include "simplnx/Parameters/GeometrySelectionParameter.hpp" @@ -76,3 +79,145 @@ TEST_CASE("SimplnxCore::ComputeTriangleGeomCentroids", "[SimplnxCore][ComputeTri UnitTest::CheckArraysInheritTupleDims(dataStructure); } + +namespace +{ +// Builds a tiny triangle mesh on a periodic domain X in [0, 4]. All vertices sit at z = 0 and each feature +// occupies a distinct y value, so ONLY the X axis can trip the periodic (both-opposing-faces) condition; +// the whole-geometry bounding box is exactly x=[0,4], y=[1,3], z=[0,0]. Face labels group triangles into +// four features that exercise the periodic centroid, plus feature 0 which is empty: +// F1 non-wrapping x = {1, 2, 3} -> touches neither x-face, centroid unchanged +// F2 symmetric wrap x = {0, 0, 4, 4} -> mass split evenly across the seam +// F3 asymmetric wrap x = {0, 3, 3.5, 4} -> bug-killer: the old constant offset lands at 4.625 (out of bounds) +// F4 domain-filling x = {0, 1, 2, 3, 4} -> degenerate, falls back to the arithmetic mean +DataStructure BuildPeriodicToyMesh() +{ + DataStructure dataStructure; + auto& triangleGeom = *TriangleGeom::Create(dataStructure, k_TriangleGeometryName); + + // 16 vertices (x, y, z). Indices 0-2 -> F1, 3-6 -> F2, 7-10 -> F3, 11-15 -> F4. + const std::vector vertices = {1.0F, 1.0F, 0.0F, 2.0F, 1.0F, 0.0F, 3.0F, 1.0F, 0.0F, // F1 + 0.0F, 2.0F, 0.0F, 0.0F, 2.0F, 0.0F, 4.0F, 2.0F, 0.0F, 4.0F, 2.0F, 0.0F, // F2 + 0.0F, 3.0F, 0.0F, 3.0F, 3.0F, 0.0F, 3.5F, 3.0F, 0.0F, 4.0F, 3.0F, 0.0F, // F3 + 0.0F, 1.0F, 0.0F, 1.0F, 1.0F, 0.0F, 2.0F, 1.0F, 0.0F, 3.0F, 1.0F, 0.0F, 4.0F, 1.0F, 0.0F}; // F4 + const usize numVertices = vertices.size() / 3; + + // 8 triangles. Each feature's triangles share vertices so the union of their vertices is the feature set. + const std::vector faces = {0, 1, 2, 3, 4, 5, 4, 5, 6, 7, 8, 9, 8, 9, 10, 11, 12, 13, 12, 13, 14, 13, 14, 15}; + const usize numFaces = faces.size() / 3; + + auto* vertexAttrMat = AttributeMatrix::Create(dataStructure, INodeGeometry0D::k_VertexAttributeMatrixName, {numVertices}, triangleGeom.getId()); + triangleGeom.setVertexAttributeMatrix(*vertexAttrMat); + auto* faceAttrMat = AttributeMatrix::Create(dataStructure, INodeGeometry2D::k_FaceAttributeMatrixName, {numFaces}, triangleGeom.getId()); + triangleGeom.setFaceAttributeMatrix(*faceAttrMat); + + auto vertexStore = std::make_unique>(std::vector{numVertices}, std::vector{3}, 0.0F); + auto* vertexList = IGeometry::SharedVertexList::Create(dataStructure, "Vertices", std::move(vertexStore), vertexAttrMat->getId()); + auto faceStore = std::make_unique>(std::vector{numFaces}, std::vector{3}, 0); + auto* faceList = IGeometry::SharedFaceList::Create(dataStructure, "SharedTriList", std::move(faceStore), faceAttrMat->getId()); + + auto vertexStoreRef = vertexList->getDataStorePtr().lock(); + for(usize i = 0; i < vertices.size(); i++) + { + vertexStoreRef->setValue(i, vertices[i]); + } + auto faceStoreRef = faceList->getDataStorePtr().lock(); + for(usize i = 0; i < faces.size(); i++) + { + faceStoreRef->setValue(i, faces[i]); + } + triangleGeom.setVertices(*vertexList); + triangleGeom.setFaceList(*faceList); + + // FaceLabels (2 components/triangle). Second label is 0 (ignored by the filter). Feature IDs 1..4. + const std::vector faceLabels = {1, 0, 2, 0, 2, 0, 3, 0, 3, 0, 4, 0, 4, 0, 4, 0}; + auto* faceDataAttrMat = AttributeMatrix::Create(dataStructure, k_FaceAttributeMatrixName, {numFaces}, triangleGeom.getId()); + auto faceLabelStore = std::make_unique>(std::vector{numFaces}, std::vector{2}, 0); + auto* faceLabelArray = Int32Array::Create(dataStructure, k_FaceLabelsName, std::move(faceLabelStore), faceDataAttrMat->getId()); + auto faceLabelStoreRef = faceLabelArray->getDataStorePtr().lock(); + for(usize i = 0; i < faceLabels.size(); i++) + { + faceLabelStoreRef->setValue(i, faceLabels[i]); + } + + // Feature AttributeMatrix sized for feature IDs 0..4 so no resize happens during execute. + AttributeMatrix::Create(dataStructure, k_FaceFeatureName, {5}, triangleGeom.getId()); + + return dataStructure; +} +} // namespace + +TEST_CASE("SimplnxCore::ComputeTriangleGeomCentroids: Periodic Minimum-Image Oracle", "[SimplnxCore][ComputeTriangleGeomCentroids]") +{ + UnitTest::LoadPlugins(); + + const DataPath centroidsPath = k_FeatureAttributeMatrixPath.createChildPath(k_CentroidsArrayName); + constexpr float32 k_Tol = 1.0e-4F; + + // Expected centroids per feature (feature 0 is empty -> default 0,0,0). + // Non-periodic: plain arithmetic mean of each feature's unique vertices. + const std::vector> nonPeriodicExpected = { + {0.0F, 0.0F, 0.0F}, // F0 empty + {2.0F, 1.0F, 0.0F}, // F1 + {2.0F, 2.0F, 0.0F}, // F2 naive x = (0+0+4+4)/4 + {2.625F, 3.0F, 0.0F}, // F3 naive x = (0+3+3.5+4)/4 + {2.0F, 1.0F, 0.0F}}; // F4 naive x = (0+1+2+3+4)/5 + + // Periodic: X component becomes the minimum-image mean on wrapping features; y/z unchanged. + const std::vector> periodicExpected = { + {0.0F, 0.0F, 0.0F}, // F0 empty + {2.0F, 1.0F, 0.0F}, // F1 does not span -> unchanged + {0.0F, 2.0F, 0.0F}, // F2 symmetric wrap -> seam at x=0 + {3.625F, 3.0F, 0.0F}, // F3 asymmetric wrap -> largest-gap mean (old code -> 4.625, out of bounds) + {2.0F, 1.0F, 0.0F}}; // F4 domain-filling -> arithmetic-mean fallback + + const bool isPeriodic = GENERATE(false, true); + const auto& expected = isPeriodic ? periodicExpected : nonPeriodicExpected; + + DYNAMIC_SECTION("IsPeriodic = " << (isPeriodic ? "true" : "false")) + { + DataStructure dataStructure = BuildPeriodicToyMesh(); + + ComputeTriangleGeomCentroidsFilter filter; + Arguments args; + args.insertOrAssign(ComputeTriangleGeomCentroidsFilter::k_TriGeometryDataPath_Key, std::make_any(k_GeometryPath)); + args.insertOrAssign(ComputeTriangleGeomCentroidsFilter::k_FaceLabelsArrayPath_Key, std::make_any(k_FaceLabelsPath)); + args.insertOrAssign(ComputeTriangleGeomCentroidsFilter::k_FeatureAttributeMatrixPath_Key, std::make_any(k_FeatureAttributeMatrixPath)); + args.insertOrAssign(ComputeTriangleGeomCentroidsFilter::k_CentroidsArrayName_Key, std::make_any(k_CentroidsArrayName)); + args.insertOrAssign(ComputeTriangleGeomCentroidsFilter::k_IsPeriodic_Key, std::make_any(isPeriodic)); + + auto preflightResult = filter.preflight(dataStructure, args); + SIMPLNX_RESULT_REQUIRE_VALID(preflightResult.outputActions); + + auto executeResult = filter.execute(dataStructure, args); + SIMPLNX_RESULT_REQUIRE_VALID(executeResult.result); + + REQUIRE_NOTHROW(dataStructure.getDataRefAs(centroidsPath)); + const auto& centroids = dataStructure.getDataRefAs(centroidsPath); + + const auto& boundingBox = dataStructure.getDataRefAs(k_GeometryPath).getBoundingBox(); + + for(usize feature = 0; feature < expected.size(); feature++) + { + for(usize comp = 0; comp < 3; comp++) + { + INFO(fmt::format("feature {} component {}", feature, comp)); + REQUIRE(std::abs(centroids[feature * 3 + comp] - expected[feature][comp]) < k_Tol); + } + } + + // Class 4 invariant: every non-empty feature centroid must lie inside the periodic domain. + for(usize feature = 1; feature < expected.size(); feature++) + { + for(usize comp = 0; comp < 3; comp++) + { + const float32 value = centroids[feature * 3 + comp]; + INFO(fmt::format("in-bounds feature {} component {}", feature, comp)); + REQUIRE(value >= boundingBox.getMinPoint()[comp] - k_Tol); + REQUIRE(value <= boundingBox.getMaxPoint()[comp] + k_Tol); + } + } + + UnitTest::CheckArraysInheritTupleDims(dataStructure); + } +} diff --git a/src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md b/src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md new file mode 100644 index 0000000000..1063807165 --- /dev/null +++ b/src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md @@ -0,0 +1,84 @@ +# V&V Report: ComputeTriangleGeomCentroidsFilter + +| | | +|--------|--------------| +| Plugin | SimplnxCore | +| SIMPLNX UUID | 074c0523-ab7a-4097-b0c3-c4dcbfb9a05e | +| DREAM3D 6.5.171 equivalent | FindTriangleGeomCentroids (DREAM3DReview) — no `FromSIMPLJson` conversion in SimplnxCore | +| Verified commit | ** | +| Status | READY FOR REVIEW | +| Sign-off | *pending second-engineer review* | + +## At a glance + +A scannable dashboard for reviewers. Each row is one sentence to one short paragraph — enough that a reader can decide whether they need to read the long-form sections below. + +| Aspect | Current state | +|------------------------|------------------------------------------------------------------------------------------------------------------------------| +| Algorithm Relationship | **Minor changes** — Port of the non-periodic core from legacy `FindTriangleGeomCentroids` (DREAM3DReview); `Is Periodic` is a SIMPLNX-only addition, now reimplemented as a minimum-image mean. | +| Oracle (confirmed) | **Class 1 + Class 4** — hand-built wrapping mesh with closed-form minimum-image centroids; encoded in `ComputeTriangleGeomCentroidsTest.cpp::"Periodic Minimum-Image Oracle"`, passes in-core **and out-of-core**. | +| Code paths enumerated | **8 of 9 exercised**; the feature-AM resize guard (path 9) is left untested (pre-sized AM, low risk). | +| Tests today | **2 test cases** — the kept `12_IN625_GBCD` exemplar (non-periodic) + the new Class 1/4 periodic oracle (`GENERATE` over `IsPeriodic`). | +| Exemplar archive | `12_IN625_GBCD.tar.gz` — non-periodic exemplar inputs+output (kept). Periodic oracle is inline (no archive). | +| Legacy comparison | **Not run yet** — legacy `FindTriangleGeomCentroids` has no periodic mode; non-periodic diff pending Step 8. | +| Bug flags | **2 SIMPLNX bugs found & fixed:** (1) constant half-domain periodic offset (issue #1665) → minimum-image mean; (2) `m_MessageHandler.m_Callback()` direct call → `std::bad_function_call` when periodic branch fired. | +| V&V phase | **All Triangle-filter V&V work complete:** discovery, Class 1/4 oracle, SIMPLNX-vs-oracle reconciliation + 2 bug fixes, algorithm review, dual-build (in-core + OOC) test pass, legacy diff (2 deviations), docs, provenance. **Status:** READY FOR REVIEW, pending second-engineer sign-off at PR. **Outstanding (separate scope):** sibling `ComputeFeatureCentroids` `ImageGeom` overload carries the same constant-offset defect (see #1665 / task 8). | + +For worked instances see `src/Plugins/OrientationAnalysis/vv/BadDataNeighborOrientationCheckFilter.md` and `src/Plugins/OrientationAnalysis/vv/ComputeAvgCAxesFilter.md` (on `topic/vv/compute_avg_caxis`). + +## Summary + +`ComputeTriangleGeomCentroids` computes the centroid of each surface-mesh feature as the mean position of the unique vertices of the triangles carrying that feature's face label, with an optional `Is Periodic` mode for features that wrap the domain boundary. V&V used a **Class 1 (Analytical)** hand-built wrapping mesh (four features: non-wrapping, symmetric wrap, asymmetric wrap, domain-filling) with closed-form minimum-image centroids, plus a **Class 4 (Invariant)** in-domain bound. The periodic path was found to use a **constant half-domain offset** that is only correct for mass symmetric about the seam (issue #1665); it was reimplemented as a largest-empty-gap minimum-image mean. A second latent bug — an unguarded `m_Callback` invocation that threw `std::bad_function_call` when the periodic branch fired — was also fixed. + +## Algorithm Relationship + +**Minor changes** (Port of the non-periodic core + SIMPLNX-only periodic addition) + +*Evidence:* No `FromSIMPLJson` conversion / not in `SimplnxCoreLegacyUUIDMapping`; the non-periodic centroid (arithmetic mean of a feature's unique vertices) mirrors legacy `FindTriangleGeomCentroids` (DREAM3DReview). The `Is Periodic` option and its wrap adjustment are SIMPLNX-only, so there is no legacy periodic behavior to preserve — the periodic path is judged on correctness alone. *(Legacy inspection pending Step 8.)* + +## Oracle + +*Class:* **1 (Analytical)** + **4 (Invariant)** + +*Applied:* A hand-built triangle mesh on periodic domain X∈[0,4] (all z=0, one distinct y per feature, so only X trips the periodic condition). Expected minimum-image centroids derived by hand via the largest-empty-gap method: F1 non-wrapping (2,1,0); F2 symmetric wrap (0,2,0); F3 asymmetric wrap (3.625,3,0) — **the discriminating case: the old constant offset yields 4.625, outside the domain**; F4 domain-filling → arithmetic-mean fallback (2,1,0). Class 4: every non-empty centroid lies within the geometry bounding box. + +*Encoded:* `test/ComputeTriangleGeomCentroidsTest.cpp::"Periodic Minimum-Image Oracle"` — `GENERATE(false, true)` over `IsPeriodic`, 4 features × 3 components + in-bounds invariant. Passes **in-core and out-of-core** (`simplnx-ooc-Rel`). + +*Second-engineer review:* *pending — the oracle design + hand-derived F3=3.625 value were confirmed with the requesting engineer before encoding; formal second-engineer review at PR.* + +## Code path coverage + +Source: `src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeTriangleGeomCentroids.cpp` (~108 lines) + the `AdjustCentroidsForPeriodicFaces` `BoundingBox` overload in `src/simplnx/Utilities/GeometryHelpers.cpp`. Phases: (a) gather per-feature vertex sets from face labels, (b) per-feature arithmetic centroid, (c) optional periodic minimum-image adjustment. + +| # | Phase | Path | Test case | +|---|-------|------|-----------| +| 1 | (a) gather | `faceLabel > 0` inserts triangle vertices into the feature set | both oracle sections (all features) | +| 2 | (a) gather | `faceLabel == 0/negative` ignored | Oracle — second label 0 on every triangle; feature 0 stays empty | +| 3 | (b) finalize | non-empty feature → arithmetic mean | `IsPeriodic=false` section (F1–F4) | +| 4 | (b) finalize | empty feature set → default (0,0,0) | Oracle — F0 asserted (0,0,0) | +| 5 | (c) periodic | axis spans both faces → minimum-image mean | `IsPeriodic=true` — F2 (unique wrap gap), F3 (unique interior gap) | +| 6 | (c) periodic | axis does not span → component unchanged | `IsPeriodic=true` — F1, and y/z of F2–F4 | +| 7 | (c) periodic | domain-filling (no dominant gap) → arithmetic fallback | `IsPeriodic=true` — F4 | +| 8 | (c) periodic | any axis adjusted → emit info message | Exercised by F2–F4 (regression pin for the `bad_function_call` fix) | +| 9 | resize | feature AM smaller than maxFeatureId+1 → resize | *Not directly tested — the exemplar test's AM is pre-sized; low-risk resize guard.* | + +## Test inventory + +| Test case | Status | Notes | +|-----------|--------|-------| +| `SimplnxCore::ComputeTriangleGeomCentroids` (exemplar, `12_IN625_GBCD`) | kept | Non-periodic exemplar comparison; unchanged. | +| `SimplnxCore::ComputeTriangleGeomCentroids: Periodic Minimum-Image Oracle` | new-for-V&V | Class 1 + Class 4; covers periodic + non-periodic; pins both bug fixes. | + +## Exemplar archive + +- **Archive:** `12_IN625_GBCD.tar.gz` (pre-existing; declared in `src/Plugins/OrientationAnalysis/test/CMakeLists.txt`, shared with several surface-mesh tests). Not regenerated by this V&V. +- **SHA512:** `f696a8af181505947e6fecfdb1a11fda6c762bba5e85fea8d484b1af00bf18643e1d930d48f092ee238d1c19c9ce7c4fb5a8092d17774bda867961a1400e9cea` +- **Provenance:** `src/Plugins/SimplnxCore/vv/provenance/12_IN625_GBCD.md` +- **Note:** The archive's `Centroids` reference array is a legacy-derived regression pin (circular oracle) for the non-periodic path only. The authoritative oracle is the inline Class 1 + Class 4 fixture; the periodic branch is not exercised by this archive. + +## Deviations from DREAM3D 6.5.171 + +Comparison method: source inspection against legacy `FindTriangleGeomCentroids` (SurfaceMeshing plugin). The non-periodic path is a line-for-line match (identical `> 0` label guard, identical float32 arithmetic-mean of set-ordered vertices), so it is bit-identical; two behavioral deviations are documented. + +- `ComputeTriangleGeomCentroidsFilter-D1` — `Is Periodic` is a SIMPLNX-only mode; its wrapped centroid (minimum-image mean) has no 6.5.171 analog — see `vv/deviations/ComputeTriangleGeomCentroidsFilter.md` +- `ComputeTriangleGeomCentroidsFilter-D2` — out-of-range face label: 6.5.171 errors (`-99500`), SIMPLNX resizes the Feature Attribute Matrix — see `vv/deviations/ComputeTriangleGeomCentroidsFilter.md` diff --git a/src/Plugins/SimplnxCore/vv/deviations/ComputeTriangleGeomCentroidsFilter.md b/src/Plugins/SimplnxCore/vv/deviations/ComputeTriangleGeomCentroidsFilter.md new file mode 100644 index 0000000000..3369ef4199 --- /dev/null +++ b/src/Plugins/SimplnxCore/vv/deviations/ComputeTriangleGeomCentroidsFilter.md @@ -0,0 +1,49 @@ +# Deviations from DREAM3D 6.5.171: ComputeTriangleGeomCentroidsFilter + +This file lists every documented behavioral difference between this SIMPLNX filter and its DREAM3D 6.5.171 equivalent. + +Entries are referenced by stable ID (`ComputeTriangleGeomCentroidsFilter-D`) from the V&V report and from public migration guidance. The ID is stable across renames; the Filter UUID field is the permanent cross-reference anchor. + +Legacy equivalent: `FindTriangleGeomCentroids` (DREAM3D 6.5.171, SurfaceMeshing plugin). Comparison method: source inspection of `FindTriangleGeomCentroids.cpp` against the SIMPLNX algorithm; the non-periodic path is a line-for-line match, so no empirical binary A/B was required (it can be run for bit-confirmation if desired). + +--- + +## ComputeTriangleGeomCentroidsFilter-D1 + +| Field | Value | +|---|---| +| **Deviation ID** | `ComputeTriangleGeomCentroidsFilter-D1` | +| **Filter UUID** | `074c0523-ab7a-4097-b0c3-c4dcbfb9a05e` | +| **Status** | active | + +**Symptom:** SIMPLNX exposes an `Is Periodic` option that DREAM3D 6.5.171 does not; when enabled, features that wrap the domain boundary receive a different (and correct) centroid than the naive arithmetic mean. + +**Root cause:** Algorithmic choice (SIMPLNX-only feature). Legacy `FindTriangleGeomCentroids` has no periodic mode — it always reports the plain arithmetic mean of a feature's unique vertices. SIMPLNX adds `Is Periodic`; for any axis on which a feature touches both opposing periodic faces, the centroid component is computed as a minimum-image (largest-empty-gap) mean of the feature's vertex coordinates (`GeometryHelpers.cpp`, `AdjustCentroidsForPeriodicFaces` `BoundingBox3Df` overload) instead of the arithmetic mean. With `Is Periodic = false` (the default) SIMPLNX reproduces the legacy result. + +**Affected users:** Only users who enable `Is Periodic` on a periodic surface mesh with grains that wrap the boundary. The default (off) path is unaffected and matches 6.5.171. + +**Recommendation:** Trust SIMPLNX. 6.5.171 has no comparable output; the SIMPLNX periodic centroid is verified against an independent Class 1 analytical oracle (see the report's Oracle section). + +--- + +## ComputeTriangleGeomCentroidsFilter-D2 + +| Field | Value | +|---|---| +| **Deviation ID** | `ComputeTriangleGeomCentroidsFilter-D2` | +| **Filter UUID** | `074c0523-ab7a-4097-b0c3-c4dcbfb9a05e` | +| **Status** | active | + +**Symptom:** When a face label is larger than the number of tuples in the target Feature Attribute Matrix, 6.5.171 aborts with error `-99500`; SIMPLNX instead grows the Feature Attribute Matrix to fit and completes successfully. + +**Root cause:** Algorithmic choice. Legacy sizes `vertexSets` to `numFeatures + 1` and errors (`-99500`) if any label meets or exceeds that bound (`FindTriangleGeomCentroids.cpp:185-192`). SIMPLNX resizes the Feature Attribute Matrix to `maxFeatureId + 1` tuples up front (`ComputeTriangleGeomCentroids.cpp:48-52`) and proceeds. The centroid values that both produce for in-range labels are identical. + +**Affected users:** Anyone who supplied a Feature Attribute Matrix smaller than `max(faceLabel) + 1`. Under 6.5.171 the run failed; under SIMPLNX it succeeds with an enlarged matrix. + +**Recommendation:** Either acceptable. SIMPLNX is strictly more permissive; the computed centroids agree for all valid labels. Users porting a pipeline that relied on the `-99500` error as a validation gate should add an explicit size check. + +--- + +## Note on the periodic bug fix (internal, not a legacy deviation) + +Prior to this V&V, SIMPLNX's periodic path added a **constant `|max − min| / 2` offset** to the naive centroid (issue #1665). That is only correct for a feature whose mass is symmetric about the seam and can place the centroid outside the domain for asymmetric wrapped features. It was corrected to the minimum-image mean described in D1. This is a SIMPLNX-internal correction, **not** a deviation from 6.5.171 (which has no periodic path), so it is recorded here only for provenance. diff --git a/src/Plugins/SimplnxCore/vv/provenance/12_IN625_GBCD.md b/src/Plugins/SimplnxCore/vv/provenance/12_IN625_GBCD.md new file mode 100644 index 0000000000..2453876024 --- /dev/null +++ b/src/Plugins/SimplnxCore/vv/provenance/12_IN625_GBCD.md @@ -0,0 +1,40 @@ +# Exemplar Archive Provenance: 12_IN625_GBCD.tar.gz + +This sidecar records how an exemplar archive used in unit tests was generated. It is the answer to "where did this gold-standard data come from?" + +--- + +## Archive identity + +| Field | Value | +|---|---| +| **Archive** | `12_IN625_GBCD.tar.gz` | +| **SHA512** | `f696a8af181505947e6fecfdb1a11fda6c762bba5e85fea8d484b1af00bf18643e1d930d48f092ee238d1c19c9ce7c4fb5a8092d17774bda867961a1400e9cea` | +| **Used by tests** | `SimplnxCore::ComputeTriangleGeomCentroids` (non-periodic exemplar); also shared by `ComputeTriangleGeomShapes`, `ComputeTriangleGeomVolumes`, `SharedFeatureFace`. | +| **Generated by** | Pre-existing BlueQuartz IN625 GBCD dataset (predates this V&V). | +| **Generated on** | Pre-existing. | +| **Generated at commit** | Pre-existing (declared in `src/Plugins/OrientationAnalysis/test/CMakeLists.txt`; not re-generated by this V&V). | + +## How it was generated + +The archive is a pre-existing IN625 grain-boundary-character-distribution surface mesh shipped in the BlueQuartz Data_Archive. This V&V did **not** regenerate or modify it. The `ComputeTriangleGeomCentroids` unit test loads it, runs the filter with `Is Periodic = false`, and compares the generated `Centroids [NX Computed]` array against the `Centroids` array already stored in the file. + +## Canonical oracle output + +| DataPath | Source of expected values | +|---|---| +| `.../FaceFeatureData/Centroids` (reference array in the archive) | **Regression pin only** — this reference array is prior DREAM3D output (a circular oracle for the non-periodic path). | + +The **primary** oracle for this filter is **not** this archive. Correctness of both the non-periodic and periodic paths is established independently by the inline **Class 1 (Analytical) + Class 4 (Invariant)** fixtures in `test/ComputeTriangleGeomCentroidsTest.cpp::"Periodic Minimum-Image Oracle"`. This exemplar test is retained only as a large-mesh regression pin for the non-periodic path; it does not exercise the periodic branch. + +## Oracle provenance (Classes 2, 3, 5 only) + +Not applicable. The independent oracle is Class 1 / Class 4 and lives directly in the test code (no archive-based oracle). + +## Second-engineer oracle review + +Pending — combined with the V&V PR review (see report Oracle section). + +## Regenerated to fix a circular-oracle situation? + +No. The archive is unchanged. The circular reference `Centroids` array is explicitly demoted to a regression pin here; the authoritative oracle is the inline analytical fixture, so no regeneration was necessary. diff --git a/src/simplnx/Utilities/GeometryHelpers.cpp b/src/simplnx/Utilities/GeometryHelpers.cpp index 62f10cd0cd..3b8f754c68 100644 --- a/src/simplnx/Utilities/GeometryHelpers.cpp +++ b/src/simplnx/Utilities/GeometryHelpers.cpp @@ -1,5 +1,11 @@ #include "GeometryHelpers.hpp" +#include +#include +#include +#include +#include + using namespace nx::core; namespace nx::core::GeometryHelpers::Description { @@ -217,25 +223,126 @@ BoundingBoxFaces FindElementPeriodicFaces(const BoundingBox3Df& boundingBox, con return edgeFaces; } -bool AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const BoundingBoxFaces& faces, Float32AbstractDataStore& centroids, IGeometry::MeshIndexType featureId) +namespace { - bool isPeriodic = false; - const auto minPoint = boundingBox.getMinPoint(); - const auto maxPoint = boundingBox.getMaxPoint(); +// Computes the periodic (minimum-image) mean of a set of 1-D coordinates on a domain of length +// domainLength whose lower bound is domainMin, using the largest-empty-gap method. Every coordinate is +// reduced into [0, domainLength); the largest circular gap between consecutive sorted values locates the +// empty region of the wrapped feature; the coordinates on the near side of that gap are unwrapped +// (+domainLength) so the feature becomes contiguous, the mean is taken, and the result is mapped back into +// [domainMin, domainMin + domainLength). When no single gap dominates (the feature fills the domain, so the +// largest gap is not unique) the plain arithmetic mean of the original coordinates is returned as a +// fallback, because the circular mean of a domain-filling distribution is not meaningful. Unlike a constant +// half-domain offset, this depends on where the feature's mass actually sits, so it is correct for +// asymmetric wrapped features and never lands outside the domain. +float32 PeriodicMean1D(const std::vector& coords, float32 domainMin, float32 domainLength) +{ + const usize numCoords = coords.size(); + if(numCoords == 0) + { + return domainMin; + } + const float32 arithmeticMean = std::accumulate(coords.cbegin(), coords.cend(), 0.0f) / static_cast(numCoords); + if(domainLength <= 0.0f) + { + return arithmeticMean; + } - if(faces.contains(BoundingBox3Df::faces_enum::left) && faces.contains(BoundingBox3Df::faces_enum::right)) + // Reduce every coordinate into [0, domainLength). A coordinate on the far face maps to 0, which is + // correct: under periodicity the two opposing faces are the same location. + std::vector reduced(numCoords); + for(usize i = 0; i < numCoords; i++) { - centroids[3 * featureId + 0] += std::abs(maxPoint[0] - minPoint[0]) / 2.0f; - isPeriodic = true; + float32 value = std::fmod(coords[i] - domainMin, domainLength); + if(value < 0.0f) + { + value += domainLength; + } + reduced[i] = value; } - if(faces.contains(BoundingBox3Df::faces_enum::top) && faces.contains(BoundingBox3Df::faces_enum::bottom)) + std::sort(reduced.begin(), reduced.end()); + + // Find the largest circular gap between consecutive sorted coordinates (the last gap wraps the ring). + float32 maxGap = -1.0f; + usize cutIndex = 0; + for(usize i = 0; i < numCoords; i++) { - centroids[3 * featureId + 1] += std::abs(maxPoint[1] - minPoint[1]) / 2.0f; - isPeriodic = true; + const float32 next = (i + 1 < numCoords) ? reduced[i + 1] : reduced[0] + domainLength; + const float32 gap = next - reduced[i]; + if(gap > maxGap) + { + maxGap = gap; + cutIndex = i; + } + } + + // If more than one gap is (within tolerance) as large as the maximum, there is no dominant empty region: + // the feature fills the domain and the wrapped mean is undefined, so fall back to the arithmetic mean. + const float32 gapTolerance = domainLength * 1.0e-6f; + usize maxGapCount = 0; + for(usize i = 0; i < numCoords; i++) + { + const float32 next = (i + 1 < numCoords) ? reduced[i + 1] : reduced[0] + domainLength; + if(std::abs((next - reduced[i]) - maxGap) <= gapTolerance) + { + maxGapCount++; + } } - if(faces.contains(BoundingBox3Df::faces_enum::front) && faces.contains(BoundingBox3Df::faces_enum::back)) + if(maxGapCount > 1) { - centroids[3 * featureId + 2] += std::abs(maxPoint[2] - minPoint[2]) / 2.0f; + return arithmeticMean; + } + + // Unwrap the coordinates on the near side of the largest gap by one domain length so the feature is + // contiguous, average in double precision, then map the result back into the domain. + float64 sum = 0.0; + for(usize i = 0; i < numCoords; i++) + { + float64 value = static_cast(reduced[i]); + if(i <= cutIndex) + { + value += static_cast(domainLength); + } + sum += value; + } + float64 mean = std::fmod(sum / static_cast(numCoords), static_cast(domainLength)); + if(mean < 0.0) + { + mean += static_cast(domainLength); + } + return domainMin + static_cast(mean); +} +} // namespace + +bool AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const BoundingBoxFaces& faces, const Float32AbstractDataStore& vertices, + const std::set& vertexSet, Float32AbstractDataStore& centroids, IGeometry::MeshIndexType featureId) +{ + bool isPeriodic = false; + const auto minPoint = boundingBox.getMinPoint(); + const auto maxPoint = boundingBox.getMaxPoint(); + + // left/right => X (axis 0), top/bottom => Y (axis 1), front/back => Z (axis 2) + const std::array, 3> axisFaces = { + std::make_pair(BoundingBox3Df::faces_enum::left, BoundingBox3Df::faces_enum::right), std::make_pair(BoundingBox3Df::faces_enum::top, BoundingBox3Df::faces_enum::bottom), + std::make_pair(BoundingBox3Df::faces_enum::front, BoundingBox3Df::faces_enum::back)}; + + // For each axis on which the feature touches both opposing periodic faces, the naive arithmetic centroid + // lands in the empty middle of the wrapped feature. Replace that component with the minimum-image mean + // computed directly from the feature's vertex coordinates on that axis. + for(usize axis = 0; axis < 3; axis++) + { + if(!faces.contains(axisFaces[axis].first) || !faces.contains(axisFaces[axis].second)) + { + continue; + } + const float32 domainLength = std::abs(maxPoint[axis] - minPoint[axis]); + std::vector axisCoords; + axisCoords.reserve(vertexSet.size()); + for(const auto& vert : vertexSet) + { + axisCoords.push_back(vertices[3 * vert + axis]); + } + centroids[3 * featureId + axis] = PeriodicMean1D(axisCoords, minPoint[axis], domainLength); isPeriodic = true; } diff --git a/src/simplnx/Utilities/GeometryHelpers.hpp b/src/simplnx/Utilities/GeometryHelpers.hpp index 49de30a025..f8799e8b2f 100644 --- a/src/simplnx/Utilities/GeometryHelpers.hpp +++ b/src/simplnx/Utilities/GeometryHelpers.hpp @@ -961,16 +961,22 @@ using BoundingBoxFaces = std::unordered_set; BoundingBoxFaces SIMPLNX_EXPORT FindElementPeriodicFaces(const BoundingBox3Df& boundingBox, const Float32AbstractDataStore& vertices, const std::set& vertexSet); /** - * @brief Adjusts centroids for periodic edge cases. The data is assumed to - * match the specified bounding box in shape as abstract shapes are not - * supported by this function. - * Returns true if the feature ID is periodic. Otherwise, returns false. - * @param boundingBox - * @param faces - * @param centroids - * @param featureId + * @brief Adjusts a feature's centroid for periodic (wrap-around) boundaries. For every axis on which the + * feature touches both opposing periodic faces, the centroid component is replaced with the minimum-image + * (largest-empty-gap) mean of the feature's vertex coordinates on that axis, rather than the naive + * arithmetic mean, which for a wrapped feature lands in the empty middle. Non-spanning axes are left + * untouched, and a domain-filling feature falls back to the arithmetic mean. The data is assumed to match + * the specified bounding box in shape as abstract shapes are not supported by this function. + * Returns true if the feature ID is periodic (any axis adjusted). Otherwise, returns false. + * @param boundingBox The bounding box of the whole geometry (the periodic domain). + * @param faces The set of bounding-box faces the feature's vertices touch. + * @param vertices The shared vertex coordinate store of the geometry. + * @param vertexSet The set of vertex indices belonging to this feature. + * @param centroids The per-feature centroid store to update (modified in place). + * @param featureId The feature whose centroid is adjusted. */ -bool SIMPLNX_EXPORT AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const BoundingBoxFaces& faces, Float32AbstractDataStore& centroids, IGeometry::MeshIndexType featureId); +bool SIMPLNX_EXPORT AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const BoundingBoxFaces& faces, const Float32AbstractDataStore& vertices, + const std::set& vertexSet, Float32AbstractDataStore& centroids, IGeometry::MeshIndexType featureId); bool SIMPLNX_EXPORT AdjustCentroidsForPeriodicFaces(const ImageGeom& imageGeom, const UInt64AbstractDataStore& xRanges, const UInt64AbstractDataStore& yRanges, const UInt64AbstractDataStore& zRanges, Float32AbstractDataStore& centroids); From 32f1aa0d4c607150a046f7c62709fd7bf6620ef6 Mon Sep 17 00:00:00 2001 From: Michael Jackson Date: Mon, 13 Jul 2026 13:20:07 -0400 Subject: [PATCH 2/2] BUG: Fix ComputeFeatureCentroids periodic centroid model (issue #1665) The periodic path applied a constant half-domain offset to the naive centroid. That model is only correct when a wrapped feature's mass is symmetric about the seam and can place asymmetric features outside the domain. The prior spacing fix (#1658) was necessary but not sufficient. * Reimplement the periodic path as a minimum-image circular mean: accumulate per-feature unit-vector (cos/sin) sums during the cell sweep, and at finalize replace the arithmetic centroid with the circular mean on axes where the feature spans the full extent. A near-zero resultant (domain-filling feature) falls back to the arithmetic mean. Trig is only evaluated when Is Periodic is enabled. * Remove the range-based ImageGeom overload of AdjustCentroidsForPeriodicFaces: per-feature index ranges are mathematically insufficient to compute a wrapped centroid, and it now has no callers. * Re-derive periodic test fixtures D and E against an independent minimum-image oracle (no longer the circular constant-offset oracle), and add Fixture F for the domain-filling fallback. Pass in-core and OOC. * Update the ComputeFeatureCentroids V&V report, deviations, and user docs to reflect the corrected model. Signed-off-by: Michael Jackson --- .../docs/ComputeFeatureCentroidsFilter.md | 2 +- .../Algorithms/ComputeFeatureCentroids.cpp | 102 +++++++++++++++++- .../test/ComputeFeatureCentroidsTest.cpp | 40 ++++--- .../test/ComputeTriangleGeomCentroidsTest.cpp | 28 +++-- .../vv/ComputeFeatureCentroidsFilter.md | 17 ++- .../vv/ComputeTriangleGeomCentroidsFilter.md | 2 +- .../ComputeFeatureCentroidsFilter.md | 4 +- src/simplnx/Utilities/GeometryHelpers.cpp | 50 +-------- src/simplnx/Utilities/GeometryHelpers.hpp | 3 - 9 files changed, 160 insertions(+), 88 deletions(-) diff --git a/src/Plugins/SimplnxCore/docs/ComputeFeatureCentroidsFilter.md b/src/Plugins/SimplnxCore/docs/ComputeFeatureCentroidsFilter.md index 17b91c122c..5dc4fd45ff 100644 --- a/src/Plugins/SimplnxCore/docs/ComputeFeatureCentroidsFilter.md +++ b/src/Plugins/SimplnxCore/docs/ComputeFeatureCentroidsFilter.md @@ -6,7 +6,7 @@ Generic (Morphological) ## Description -This **Filter** calculates the *centroid* of each **Feature** by determining the average X, Y, and Z position (in physical coordinates) of all the **Cells** belonging to the **Feature**. The per-cell coordinates are accumulated using Kahan compensated summation to limit floating-point round-off on features with large cell counts. An *Is Periodic* option is available: when enabled, a **Feature** that spans the full extent of an axis (touching both opposing faces of the **Image Geometry**) is treated as wrapping around to the opposite face, and its centroid on that axis is shifted by half the physical distance between the first and last cell centers. When *Is Periodic* is disabled, **Features** that intersect the outer surfaces of the sample will still have *centroids* calculated, but they will be *centroids* of the truncated part of the **Feature** that lies inside the sample. +This **Filter** calculates the *centroid* of each **Feature** by determining the average X, Y, and Z position (in physical coordinates) of all the **Cells** belonging to the **Feature**. The per-cell coordinates are accumulated using Kahan compensated summation to limit floating-point round-off on features with large cell counts. An *Is Periodic* option is available: when enabled, a **Feature** that spans the full extent of an axis (touching both opposing faces of the **Image Geometry**) is treated as wrapping around to the opposite face, and its centroid on that axis is computed as a minimum-image (circular) mean — the cell positions are unwrapped across the boundary before averaging and the result is mapped back into the domain. This depends on where the **Feature**'s cells actually lie, so it is correct for features whose mass wraps asymmetrically and never places the centroid outside the sample; a feature that fills the whole domain along an axis falls back to the ordinary average on that axis. When *Is Periodic* is disabled, **Features** that intersect the outer surfaces of the sample will still have *centroids* calculated, but they will be *centroids* of the truncated part of the **Feature** that lies inside the sample. A **Feature** with no **Cells** (an unused Feature Id) keeps a centroid of (0, 0, 0). diff --git a/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeFeatureCentroids.cpp b/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeFeatureCentroids.cpp index 8cbc33fe10..09c6205bbf 100644 --- a/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeFeatureCentroids.cpp +++ b/src/Plugins/SimplnxCore/src/SimplnxCore/Filters/Algorithms/ComputeFeatureCentroids.cpp @@ -4,10 +4,21 @@ #include "simplnx/DataStructure/DataGroup.hpp" #include "simplnx/DataStructure/Geometry/ImageGeom.hpp" #include "simplnx/Utilities/DataArrayUtilities.hpp" -#include "simplnx/Utilities/GeometryHelpers.hpp" #include "simplnx/Utilities/ParallelDataAlgorithm.hpp" #include +#include +#include + +namespace +{ +// 2*pi, used to map cell-center positions onto the periodic unit circle. +constexpr double k_TwoPi = 6.283185307179586; +// A feature whose per-axis unit-vector resultant falls below this length has its mass spread ~uniformly +// around the domain (a domain-filling feature); its circular mean is not meaningful, so the arithmetic +// mean is kept instead. +constexpr double k_DegenerateResultant = 1.0e-6; +} // namespace using namespace nx::core; @@ -21,7 +32,8 @@ class ComputeFeatureCentroidsImpl1 { public: ComputeFeatureCentroidsImpl1(Float64AbstractDataStore& sum, Float64AbstractDataStore& compensation, UInt64AbstractDataStore& count, std::array dims, const nx::core::ImageGeom& imageGeom, - const Int32AbstractDataStore& featureIds, UInt64AbstractDataStore& rangeXStoreRef, UInt64AbstractDataStore& rangeYStoreRef, UInt64AbstractDataStore& rangeZStoreRef) + const Int32AbstractDataStore& featureIds, UInt64AbstractDataStore& rangeXStoreRef, UInt64AbstractDataStore& rangeYStoreRef, UInt64AbstractDataStore& rangeZStoreRef, + bool isPeriodic, Float64AbstractDataStore& sumCos, Float64AbstractDataStore& sumSin, std::array origin, std::array domainLength) : m_Sum(sum) , m_Compensation(compensation) , m_Count(count) @@ -31,6 +43,11 @@ class ComputeFeatureCentroidsImpl1 , m_RangeXStoreRef(rangeXStoreRef) , m_RangeYStoreRef(rangeYStoreRef) , m_RangeZStoreRef(rangeZStoreRef) + , m_IsPeriodic(isPeriodic) + , m_SumCos(sumCos) + , m_SumSin(sumSin) + , m_Origin(origin) + , m_DomainLength(domainLength) { } ~ComputeFeatureCentroidsImpl1() = default; @@ -85,6 +102,21 @@ class ComputeFeatureCentroidsImpl1 m_Compensation[featureId_idx] = (temp - m_Sum[featureId_idx]) - componentValue; m_Sum[featureId_idx] = temp; m_Count[featureId_idx].inc(); + + // For periodic runs, also accumulate the unit vector of each cell center's angular position + // around the domain on each axis. The per-feature (cos, sin) sums produce a minimum-image + // centroid at finalize for features that wrap the boundary. + if(m_IsPeriodic) + { + const nx::core::Point3Dd center = {voxel_center[0], voxel_center[1], voxel_center[2]}; + for(size_t axis = 0; axis < 3; axis++) + { + const double phase = k_TwoPi * (center[axis] - m_Origin[axis]) / m_DomainLength[axis]; + const size_t axisIdx = featureId * 3ULL + axis; + m_SumCos[axisIdx] = m_SumCos.getValue(axisIdx) + std::cos(phase); + m_SumSin[axisIdx] = m_SumSin.getValue(axisIdx) + std::sin(phase); + } + } } } } @@ -105,6 +137,11 @@ class ComputeFeatureCentroidsImpl1 UInt64AbstractDataStore& m_RangeXStoreRef; UInt64AbstractDataStore& m_RangeYStoreRef; UInt64AbstractDataStore& m_RangeZStoreRef; + bool m_IsPeriodic = false; + Float64AbstractDataStore& m_SumCos; + Float64AbstractDataStore& m_SumSin; + std::array m_Origin = {0.0, 0.0, 0.0}; + std::array m_DomainLength = {0.0, 0.0, 0.0}; }; } // namespace @@ -169,6 +206,21 @@ Result<> ComputeFeatureCentroids::operator()() compensation.fill(0.0); count.fill(0.0); + // Per-feature/per-axis unit-vector sums for the periodic (circular-mean) centroid. Only populated when + // Is Periodic is enabled; harmless (all zero) otherwise. + auto sumCosPtr = DataStoreUtilities::CreateDataStore(tupleShape, componentShape, IDataAction::Mode::Execute); + auto sumSinPtr = DataStoreUtilities::CreateDataStore(tupleShape, componentShape, IDataAction::Mode::Execute); + Float64AbstractDataStore& sumCos = *sumCosPtr.get(); + Float64AbstractDataStore& sumSin = *sumSinPtr.get(); + sumCos.fill(0.0); + sumSin.fill(0.0); + + const auto geomOrigin = imageGeom.getOrigin(); + const auto geomSpacing = imageGeom.getSpacing(); + const std::array origin = {static_cast(geomOrigin[0]), static_cast(geomOrigin[1]), static_cast(geomOrigin[2])}; + // Periodic domain length on each axis = number of cells * spacing (the full physical extent). + const std::array domainLength = {static_cast(xPoints) * geomSpacing[0], static_cast(yPoints) * geomSpacing[1], static_cast(zPoints) * geomSpacing[2]}; + // Create data stores to check if feature IDs are periodic componentShape[0] = 2; auto rangeXStorePtr = DataStoreUtilities::CreateDataStore(tupleShape, componentShape, IDataAction::Mode::Execute); @@ -187,7 +239,8 @@ Result<> ComputeFeatureCentroids::operator()() // by the total number of cores/threads and do a ParallelTask Algorithm instead // we might see some speedup. dataAlg.setParallelizationEnabled(false); - dataAlg.execute(ComputeFeatureCentroidsImpl1(sum, compensation, count, {xPoints, yPoints, zPoints}, imageGeom, featureIdsStoreRef, rangeXStoreRef, rangeYStoreRef, rangeZStoreRef)); + dataAlg.execute(ComputeFeatureCentroidsImpl1(sum, compensation, count, {xPoints, yPoints, zPoints}, imageGeom, featureIdsStoreRef, rangeXStoreRef, rangeYStoreRef, rangeZStoreRef, + m_InputValues->IsPeriodic, sumCos, sumSin, origin, domainLength)); // Here we are only looping over the number of features so let this just go in serial mode. // The count store carries the same voxel count in all three components of a feature; a feature with @@ -216,9 +269,48 @@ Result<> ComputeFeatureCentroids::operator()() if(m_InputValues->IsPeriodic) { m_MessageHandler({IFilter::Message::Type::Info, "Checking for periodic data."}); - if(GeometryHelpers::Topology::AdjustCentroidsForPeriodicFaces(imageGeom, rangeXStoreRef, rangeYStoreRef, rangeZStoreRef, centroids)) + + // For each axis on which a feature spans the full extent (a cell at index 0 and at the last index), the + // naive arithmetic centroid lands in the empty middle of the wrapped feature. Replace that component with + // the circular (minimum-image) mean derived from the per-feature unit-vector sums. A feature whose mass + // is spread ~uniformly around the domain (near-zero resultant) keeps its arithmetic mean. + const std::array rangeStores = {&rangeXStoreRef, &rangeYStoreRef, &rangeZStoreRef}; + const std::array dims = {xPoints, yPoints, zPoints}; + bool anyAdjusted = false; + for(size_t featureId = 0; featureId < totalFeatures; featureId++) + { + for(size_t axis = 0; axis < 3; axis++) + { + const size_t axisIdx = featureId * 3 + axis; + if(count[axisIdx] == 0) + { + continue; + } + const UInt64AbstractDataStore& rangeStore = *rangeStores[axis]; + const bool spansExtent = (rangeStore.getValue(featureId * 2 + 0) == 0 && rangeStore.getValue(featureId * 2 + 1) == dims[axis] - 1); + if(!spansExtent) + { + continue; + } + const double sumCosValue = sumCos.getValue(axisIdx); + const double sumSinValue = sumSin.getValue(axisIdx); + const double resultant = std::sqrt(sumCosValue * sumCosValue + sumSinValue * sumSinValue) / static_cast(count[axisIdx]); + if(resultant < k_DegenerateResultant) + { + continue; + } + double phase = std::atan2(sumSinValue, sumCosValue); + if(phase < 0.0) + { + phase += k_TwoPi; + } + centroids[axisIdx] = static_cast(origin[axis] + (phase / k_TwoPi) * domainLength[axis]); + anyAdjusted = true; + } + } + if(anyAdjusted) { - m_MessageHandler({IFilter::Message::Type::Info, "ComputeFeatureCentroids found Non-Contiguous Features. Centroids may require additional checks."}); + m_MessageHandler({IFilter::Message::Type::Info, "ComputeFeatureCentroids adjusted centroids of features that wrap the periodic boundary."}); } } diff --git a/src/Plugins/SimplnxCore/test/ComputeFeatureCentroidsTest.cpp b/src/Plugins/SimplnxCore/test/ComputeFeatureCentroidsTest.cpp index 699607d9a1..498412ca4c 100644 --- a/src/Plugins/SimplnxCore/test/ComputeFeatureCentroidsTest.cpp +++ b/src/Plugins/SimplnxCore/test/ComputeFeatureCentroidsTest.cpp @@ -160,10 +160,10 @@ TEST_CASE("SimplnxCore::ComputeFeatureCentroidsFilter: Class 1/4 - Periodic Boun using namespace CentroidToy; UnitTest::LoadPlugins(); - SECTION("Fixture D - periodic wrap, unit spacing") + SECTION("Fixture D - periodic wrap, unit spacing (asymmetric)") { // 4x1x1; FeatureIds [1,2,1,1]; feature 1 spans x=0 and x=3 (full extent), feature 2 (x=1) does not. - // feature 1 cells x=0,2,3 -> centers 0.5,2.5,3.5 -> mean 6.5/3 = 2.16667 + // feature 1 cells x=0,2,3 -> centers 0.5,2.5,3.5 -> non-periodic mean 6.5/3 = 2.16667 auto sNP = Build(4, 1, 1, {1.0f, 1.0f, 1.0f}, {0.0f, 0.0f, 0.0f}, 3, {1, 2, 1, 1}); auto nonPeriodic = Run(sNP, false); RequireCentroid(nonPeriodic, 1, 6.5f / 3.0f, 0.5f, 0.5f); @@ -171,26 +171,38 @@ TEST_CASE("SimplnxCore::ComputeFeatureCentroidsFilter: Class 1/4 - Periodic Boun auto sP = Build(4, 1, 1, {1.0f, 1.0f, 1.0f}, {0.0f, 0.0f, 0.0f}, 3, {1, 2, 1, 1}); auto periodic = Run(sP, true); - // Class 1 (given the DREAM3D periodic-shift model): spanning feature 1 gets +(dim-1)/2 = +1.5. - RequireCentroid(periodic, 1, 6.5f / 3.0f + 1.5f, 0.5f, 0.5f); + // Class 1 minimum-image oracle: domain L=4 (seam at x=0). The cell at center 0.5 wraps to 4.5, so the + // contiguous cluster is {2.5, 3.5, 4.5}; mean = 10.5/3 = 3.5, mapped back into [0,4) -> 3.5. + // (The old constant-offset model gave 6.5/3 + 1.5 = 3.6667, which is wrong.) + RequireCentroid(periodic, 1, 3.5f, 0.5f, 0.5f); // Class 4 invariant: the non-spanning feature 2 must be unchanged by the periodic pass. RequireCentroid(periodic, 2, 1.5f, 0.5f, 0.5f); } - SECTION("Fixture E - periodic wrap, NON-unit spacing (regression pin for the spacing fix)") + SECTION("Fixture E - periodic wrap, NON-unit spacing (asymmetric, interior result)") { - // 4x1x1, spacing (2,1,1), origin (10,0,0); FeatureIds [1,2,2,1] - // x-centers: idx0 -> 11.0, idx3 -> 17.0; non-periodic centroid.x = 14.0 - auto sNP = Build(4, 1, 1, {2.0f, 1.0f, 1.0f}, {10.0f, 0.0f, 0.0f}, 3, {1, 2, 2, 1}); + // 4x1x1, spacing (2,1,1), origin (10,0,0); FeatureIds [1,1,2,1] + // feature 1 cells idx 0,1,3 -> x-centers 11.0, 13.0, 17.0; non-periodic centroid.x = 41/3 = 13.6667 + auto sNP = Build(4, 1, 1, {2.0f, 1.0f, 1.0f}, {10.0f, 0.0f, 0.0f}, 3, {1, 1, 2, 1}); auto nonPeriodic = Run(sNP, false); - RequireCentroid(nonPeriodic, 1, 14.0f, 0.5f, 0.5f); // unambiguous analytical value + RequireCentroid(nonPeriodic, 1, 41.0f / 3.0f, 0.5f, 0.5f); - auto sP = Build(4, 1, 1, {2.0f, 1.0f, 1.0f}, {10.0f, 0.0f, 0.0f}, 3, {1, 2, 2, 1}); + auto sP = Build(4, 1, 1, {2.0f, 1.0f, 1.0f}, {10.0f, 0.0f, 0.0f}, 3, {1, 1, 2, 1}); auto periodic = Run(sP, true); - // The periodic offset must scale with spacing: (dim-1)*spacing_x/2 = 3*2/2 = 3.0 -> 14.0 + 3.0 = 17.0. - // Before the fix this returned 15.5 (offset 1.5 in cell units, ignoring spacing). Regression pin for - // the AdjustCentroidsForPeriodicFaces spacing fix. See vv/ComputeFeatureCentroidsFilter.md Phase 6. - RequireCentroid(periodic, 1, 17.0f, 0.5f, 0.5f); + // Class 1 minimum-image oracle: domain L=8 (seam at x=10 == x=18). Positions relative to origin are + // {1, 3, 7}; the largest empty gap is 3->7, so 1 and 3 wrap to 9 and 11. Cluster {7, 9, 11}; mean = 9, + // mapped back into [10,18) -> x = 11.0. (The old constant-offset model gave 13.6667 + 3.0 = 16.6667.) + RequireCentroid(periodic, 1, 11.0f, 0.5f, 0.5f); + } + + SECTION("Fixture F - domain-filling feature (arithmetic-mean fallback)") + { + // 4x1x1, unit spacing; a single feature fills every cell (centers 0.5,1.5,2.5,3.5). It spans the full + // extent, so the periodic branch fires, but the unit vectors cancel (zero resultant): the circular mean + // is undefined, so the arithmetic mean 2.0 is kept. Class 4 fallback path. + auto sP = Build(4, 1, 1, {1.0f, 1.0f, 1.0f}, {0.0f, 0.0f, 0.0f}, 2, {1, 1, 1, 1}); + auto periodic = Run(sP, true); + RequireCentroid(periodic, 1, 2.0f, 0.5f, 0.5f); } } diff --git a/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp b/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp index 0258a4643e..5dacf0517c 100644 --- a/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp +++ b/src/Plugins/SimplnxCore/test/ComputeTriangleGeomCentroidsTest.cpp @@ -96,9 +96,9 @@ DataStructure BuildPeriodicToyMesh() auto& triangleGeom = *TriangleGeom::Create(dataStructure, k_TriangleGeometryName); // 16 vertices (x, y, z). Indices 0-2 -> F1, 3-6 -> F2, 7-10 -> F3, 11-15 -> F4. - const std::vector vertices = {1.0F, 1.0F, 0.0F, 2.0F, 1.0F, 0.0F, 3.0F, 1.0F, 0.0F, // F1 - 0.0F, 2.0F, 0.0F, 0.0F, 2.0F, 0.0F, 4.0F, 2.0F, 0.0F, 4.0F, 2.0F, 0.0F, // F2 - 0.0F, 3.0F, 0.0F, 3.0F, 3.0F, 0.0F, 3.5F, 3.0F, 0.0F, 4.0F, 3.0F, 0.0F, // F3 + const std::vector vertices = {1.0F, 1.0F, 0.0F, 2.0F, 1.0F, 0.0F, 3.0F, 1.0F, 0.0F, // F1 + 0.0F, 2.0F, 0.0F, 0.0F, 2.0F, 0.0F, 4.0F, 2.0F, 0.0F, 4.0F, 2.0F, 0.0F, // F2 + 0.0F, 3.0F, 0.0F, 3.0F, 3.0F, 0.0F, 3.5F, 3.0F, 0.0F, 4.0F, 3.0F, 0.0F, // F3 0.0F, 1.0F, 0.0F, 1.0F, 1.0F, 0.0F, 2.0F, 1.0F, 0.0F, 3.0F, 1.0F, 0.0F, 4.0F, 1.0F, 0.0F}; // F4 const usize numVertices = vertices.size() / 3; @@ -156,20 +156,18 @@ TEST_CASE("SimplnxCore::ComputeTriangleGeomCentroids: Periodic Minimum-Image Ora // Expected centroids per feature (feature 0 is empty -> default 0,0,0). // Non-periodic: plain arithmetic mean of each feature's unique vertices. - const std::vector> nonPeriodicExpected = { - {0.0F, 0.0F, 0.0F}, // F0 empty - {2.0F, 1.0F, 0.0F}, // F1 - {2.0F, 2.0F, 0.0F}, // F2 naive x = (0+0+4+4)/4 - {2.625F, 3.0F, 0.0F}, // F3 naive x = (0+3+3.5+4)/4 - {2.0F, 1.0F, 0.0F}}; // F4 naive x = (0+1+2+3+4)/5 + const std::vector> nonPeriodicExpected = {{0.0F, 0.0F, 0.0F}, // F0 empty + {2.0F, 1.0F, 0.0F}, // F1 + {2.0F, 2.0F, 0.0F}, // F2 naive x = (0+0+4+4)/4 + {2.625F, 3.0F, 0.0F}, // F3 naive x = (0+3+3.5+4)/4 + {2.0F, 1.0F, 0.0F}}; // F4 naive x = (0+1+2+3+4)/5 // Periodic: X component becomes the minimum-image mean on wrapping features; y/z unchanged. - const std::vector> periodicExpected = { - {0.0F, 0.0F, 0.0F}, // F0 empty - {2.0F, 1.0F, 0.0F}, // F1 does not span -> unchanged - {0.0F, 2.0F, 0.0F}, // F2 symmetric wrap -> seam at x=0 - {3.625F, 3.0F, 0.0F}, // F3 asymmetric wrap -> largest-gap mean (old code -> 4.625, out of bounds) - {2.0F, 1.0F, 0.0F}}; // F4 domain-filling -> arithmetic-mean fallback + const std::vector> periodicExpected = {{0.0F, 0.0F, 0.0F}, // F0 empty + {2.0F, 1.0F, 0.0F}, // F1 does not span -> unchanged + {0.0F, 2.0F, 0.0F}, // F2 symmetric wrap -> seam at x=0 + {3.625F, 3.0F, 0.0F}, // F3 asymmetric wrap -> largest-gap mean (old code -> 4.625, out of bounds) + {2.0F, 1.0F, 0.0F}}; // F4 domain-filling -> arithmetic-mean fallback const bool isPeriodic = GENERATE(false, true); const auto& expected = isPeriodic ? periodicExpected : nonPeriodicExpected; diff --git a/src/Plugins/SimplnxCore/vv/ComputeFeatureCentroidsFilter.md b/src/Plugins/SimplnxCore/vv/ComputeFeatureCentroidsFilter.md index e181e122ef..aa0c68453b 100644 --- a/src/Plugins/SimplnxCore/vv/ComputeFeatureCentroidsFilter.md +++ b/src/Plugins/SimplnxCore/vv/ComputeFeatureCentroidsFilter.md @@ -35,12 +35,25 @@ This is the engineer's working V&V doc for `ComputeFeatureCentroidsFilter`, prod | Tests today | **4 TEST_CASEs (5 ctest entries), all pass**: Class 1 Analytical (Fixtures A/B/C, 3 SECTIONs), Class 1/4 Periodic (D/E), validation-error (-5351), SIMPL 6.4/6.5 backwards-compat (`DYNAMIC_SECTION`). Retired the circular consistency test. | | Exemplar archive | **None — inline analytical fixtures** (provenance sidecar written). Shared `6_6_stats_test_v2.tar.gz` no longer consumed by this test (kept for 5 other tests). `6_6_find_feature_centroids.tar.gz` **kept** (used by ExtractComponentAsArray / WriteAbaqusHexahedron — not an orphan). | | Legacy comparison | **Source-inspection three-way** (6.5.171 / 6.5.172 / SIMPLNX). Non-periodic path is an exact Kahan port; 2 deviations documented (D1 precision float32→float64 voxel center; D2 periodic feature addition). Empirical binary A/B available if bit-confirmation needed. | -| Bug flags | **1 SIMPLNX bug found & FIXED (Phase 6):** `AdjustCentroidsForPeriodicFaces` (ImageGeom overload) added `(dim−1)/2` in cell-index units to a physical-coordinate centroid, ignoring spacing. Fixed to scale by spacing (`GeometryHelpers.cpp:245`); Fixture E is the regression pin. Sibling BoundingBox overload audited — sound (physical units), left unchanged. Retroactive report's "D1: Kahan-vs-naive" is **retracted** — legacy already uses Kahan. | +| Bug flags | **2 SIMPLNX bugs found & FIXED.** (Phase 6) The spacing sub-bug: `AdjustCentroidsForPeriodicFaces` added `(dim−1)/2` in cell-index units to a physical-coordinate centroid, ignoring spacing. **(Follow-up, 2026-07, issue #1665)** The spacing fix was necessary but NOT sufficient — the underlying **constant half-domain offset model is itself wrong** (correct only for mass symmetric about the seam; can place the centroid outside the domain for asymmetric wrapped features). The periodic path was reimplemented as a **minimum-image circular mean** accumulated during the cell sweep; the range-based `ImageGeom` overload of `AdjustCentroidsForPeriodicFaces` was **removed**. See the Follow-up section below. | | V&V phase | **Phases 1–11 complete.** Oracle verified, 1 bug fixed + pinned, algorithm review applied, 8/8 code paths covered, deviations + provenance written, docs updated. Outstanding: **OOC dual-build deferred by developer**; Phase 12 (OneDrive archive); Phase 13 (status promotion + retroactive-report/INDEX update); second-engineer oracle review. | +## Follow-up (2026-07, issue #1665) — periodic model reimplemented + +The original V&V (#1658) fixed a **spacing** sub-bug in the periodic path but validated the periodic result against the **constant half-domain offset model itself** (Fixture D's expected value was "naive + (dim−1)/2", a circular oracle — the test was designed around the code). An adversarial review (#1665) showed that model is fundamentally wrong: a constant offset is correct only when a wrapped feature's mass is symmetric about the seam, and for asymmetric features it places the centroid in the empty middle of the feature or outside the domain entirely. + +This follow-up (on top of merged #1658) corrects it: + +- **Periodic path reimplemented as a minimum-image circular mean.** During the existing cell sweep, per feature and axis, the unit vectors of each cell center's angular position around the domain are accumulated (`sumCos`, `sumSin`). At finalize, for any axis on which the feature spans the full extent, the centroid component becomes `atan2(sumSin, sumCos)` mapped back into the domain. A near-zero resultant (domain-filling feature) falls back to the arithmetic mean. O(numFeatures) extra memory; trig only when `Is Periodic` is on. +- **The range-based `ImageGeom` overload of `AdjustCentroidsForPeriodicFaces` was removed** — it received only per-feature index ranges, which are mathematically insufficient to compute a wrapped centroid (the whole point of #1665). Nothing else called it. +- **Fixtures D and E re-derived against an independent minimum-image oracle** (no longer circular): D → x = 3.5 (was 3.667), E (now asymmetric, non-unit spacing) → x = 11.0 (was 17.0). New Fixture F pins the domain-filling arithmetic-mean fallback. All pass in-core and out-of-core. +- The earlier "sibling BoundingBox overload audited — sound, left unchanged" conclusion is **superseded**: that overload carried the same constant-offset defect and was fixed under issue #1665 in the same branch as this follow-up. + +The non-periodic path is unchanged (still the verified Kahan port). Deviation `ComputeFeatureCentroidsFilter-D2` is updated accordingly. + ## Summary -`ComputeFeatureCentroidsFilter` computes the centroid of each feature as the average X/Y/Z position of all cells belonging to that feature, using Kahan compensated summation for numerical precision. An optional `Is Periodic` mode wraps features that span the image boundary. Verification approach (planned): **Class 1 (Analytical) toy fixtures** — the centroid is the arithmetic mean of voxel-center coordinates, hand-derivable — paired with **Class 4 (Invariant)** bounding-box checks, then a legacy diff against DREAM3D 6.5.171 that is expected to reduce to a single float32-vs-float64 precision non-deviation plus the SIMPLNX-only periodic feature. +`ComputeFeatureCentroidsFilter` computes the centroid of each feature as the average X/Y/Z position of all cells belonging to that feature, using Kahan compensated summation for numerical precision. An optional `Is Periodic` mode wraps features that span the image boundary using a minimum-image circular mean (see the Follow-up section for the 2026-07 reimplementation). Verification approach (planned): **Class 1 (Analytical) toy fixtures** — the centroid is the arithmetic mean of voxel-center coordinates, hand-derivable — paired with **Class 4 (Invariant)** bounding-box checks, then a legacy diff against DREAM3D 6.5.171 that is expected to reduce to a single float32-vs-float64 precision non-deviation plus the SIMPLNX-only periodic feature. --- diff --git a/src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md b/src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md index 1063807165..5bc26f24ce 100644 --- a/src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md +++ b/src/Plugins/SimplnxCore/vv/ComputeTriangleGeomCentroidsFilter.md @@ -22,7 +22,7 @@ A scannable dashboard for reviewers. Each row is one sentence to one short parag | Exemplar archive | `12_IN625_GBCD.tar.gz` — non-periodic exemplar inputs+output (kept). Periodic oracle is inline (no archive). | | Legacy comparison | **Not run yet** — legacy `FindTriangleGeomCentroids` has no periodic mode; non-periodic diff pending Step 8. | | Bug flags | **2 SIMPLNX bugs found & fixed:** (1) constant half-domain periodic offset (issue #1665) → minimum-image mean; (2) `m_MessageHandler.m_Callback()` direct call → `std::bad_function_call` when periodic branch fired. | -| V&V phase | **All Triangle-filter V&V work complete:** discovery, Class 1/4 oracle, SIMPLNX-vs-oracle reconciliation + 2 bug fixes, algorithm review, dual-build (in-core + OOC) test pass, legacy diff (2 deviations), docs, provenance. **Status:** READY FOR REVIEW, pending second-engineer sign-off at PR. **Outstanding (separate scope):** sibling `ComputeFeatureCentroids` `ImageGeom` overload carries the same constant-offset defect (see #1665 / task 8). | +| V&V phase | **All Triangle-filter V&V work complete:** discovery, Class 1/4 oracle, SIMPLNX-vs-oracle reconciliation + 2 bug fixes, algorithm review, dual-build (in-core + OOC) test pass, legacy diff (2 deviations), docs, provenance. **Status:** READY FOR REVIEW, pending second-engineer sign-off at PR. The sibling `ComputeFeatureCentroids` `ImageGeom` overload carried the same constant-offset defect; it was reimplemented as a circular mean and the overload removed in the same branch (see `vv/ComputeFeatureCentroidsFilter.md` Follow-up). | For worked instances see `src/Plugins/OrientationAnalysis/vv/BadDataNeighborOrientationCheckFilter.md` and `src/Plugins/OrientationAnalysis/vv/ComputeAvgCAxesFilter.md` (on `topic/vv/compute_avg_caxis`). diff --git a/src/Plugins/SimplnxCore/vv/deviations/ComputeFeatureCentroidsFilter.md b/src/Plugins/SimplnxCore/vv/deviations/ComputeFeatureCentroidsFilter.md index 3d206025c8..780c554175 100644 --- a/src/Plugins/SimplnxCore/vv/deviations/ComputeFeatureCentroidsFilter.md +++ b/src/Plugins/SimplnxCore/vv/deviations/ComputeFeatureCentroidsFilter.md @@ -36,13 +36,13 @@ Entries are referenced by stable ID (`ComputeFeatureCentroidsFilter-D`) from **Symptom:** SIMPLNX exposes an `Is Periodic` option that, when enabled, offsets the centroid of any feature spanning the full extent of an axis (wrap-aware centroid). DREAM3D 6.5.171 `FindFeatureCentroids` has no such option and always produces the literal (truncated) centroid. -**Root cause:** *Algorithmic choice (feature addition).* The `Is Periodic` parameter and the call to `GeometryHelpers::Topology::AdjustCentroidsForPeriodicFaces` are SIMPLNX-only. Legacy has no periodic path, so there is nothing to diff against — this is an additive capability, not a divergence on shared behavior. The SIMPLNX default is `false`, which reproduces legacy exactly. SIMPL 6.4/6.5 pipelines convert with `Is Periodic` defaulted to `false`, so auto-converted legacy pipelines match 6.5.171. +**Root cause:** *Algorithmic choice (feature addition).* The `Is Periodic` parameter and the periodic centroid computation are SIMPLNX-only. Legacy has no periodic path, so there is nothing to diff against — this is an additive capability, not a divergence on shared behavior. The SIMPLNX default is `false`, which reproduces legacy exactly. SIMPL 6.4/6.5 pipelines convert with `Is Periodic` defaulted to `false`, so auto-converted legacy pipelines match 6.5.171. **Affected users:** Only users who explicitly enable `Is Periodic`. Auto-converted legacy pipelines are unaffected (default `false`). **Recommendation:** *Trust SIMPLNX.* Opt-in feature with a legacy-preserving default. -**Note (SIMPLNX-internal bug fixed during V&V — not a legacy deviation):** Prior to this V&V cycle the periodic offset in the `ImageGeom` overload of `AdjustCentroidsForPeriodicFaces` was added in **cell-index units** (`(dim−1)/2`) directly to a physical-coordinate centroid, so it was only correct when spacing == 1. This was fixed to scale by spacing (`src/simplnx/Utilities/GeometryHelpers.cpp`), with `test/ComputeFeatureCentroidsTest.cpp` Fixture E (spacing = 2) as the regression pin. Because 6.5.171 has no periodic path, this fix is not a change relative to legacy; it is recorded here for completeness and for any SIMPLNX user who ran `Is Periodic` on non-unit-spacing data before the fix. +**Note (SIMPLNX-internal bugs fixed during V&V — not legacy deviations):** The periodic centroid had two SIMPLNX-internal defects, both fixed. (1) *Spacing:* the offset in the `ImageGeom` overload of `AdjustCentroidsForPeriodicFaces` was added in cell-index units (`(dim−1)/2`) to a physical-coordinate centroid, correct only when spacing == 1 (fixed during #1658). (2) *Model (issue #1665, 2026-07 follow-up):* the constant half-domain offset is itself only correct for a feature whose mass is symmetric about the seam, and can place asymmetric wrapped features outside the domain. The periodic path was reimplemented as a **minimum-image circular mean** accumulated during the cell sweep, and the range-based `ImageGeom` overload (insufficient inputs to compute a wrapped centroid) was removed. `test/ComputeFeatureCentroidsTest.cpp` Fixtures D/E (independent minimum-image oracle) and F (domain-filling fallback) pin the corrected behavior. Because 6.5.171 has no periodic path, neither fix is a change relative to legacy; recorded here for any SIMPLNX user who ran `Is Periodic` before the fixes. --- diff --git a/src/simplnx/Utilities/GeometryHelpers.cpp b/src/simplnx/Utilities/GeometryHelpers.cpp index 3b8f754c68..9a6caf57e6 100644 --- a/src/simplnx/Utilities/GeometryHelpers.cpp +++ b/src/simplnx/Utilities/GeometryHelpers.cpp @@ -314,17 +314,17 @@ float32 PeriodicMean1D(const std::vector& coords, float32 domainMin, fl } } // namespace -bool AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const BoundingBoxFaces& faces, const Float32AbstractDataStore& vertices, - const std::set& vertexSet, Float32AbstractDataStore& centroids, IGeometry::MeshIndexType featureId) +bool AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const BoundingBoxFaces& faces, const Float32AbstractDataStore& vertices, const std::set& vertexSet, + Float32AbstractDataStore& centroids, IGeometry::MeshIndexType featureId) { bool isPeriodic = false; const auto minPoint = boundingBox.getMinPoint(); const auto maxPoint = boundingBox.getMaxPoint(); // left/right => X (axis 0), top/bottom => Y (axis 1), front/back => Z (axis 2) - const std::array, 3> axisFaces = { - std::make_pair(BoundingBox3Df::faces_enum::left, BoundingBox3Df::faces_enum::right), std::make_pair(BoundingBox3Df::faces_enum::top, BoundingBox3Df::faces_enum::bottom), - std::make_pair(BoundingBox3Df::faces_enum::front, BoundingBox3Df::faces_enum::back)}; + const std::array, 3> axisFaces = {std::make_pair(BoundingBox3Df::faces_enum::left, BoundingBox3Df::faces_enum::right), + std::make_pair(BoundingBox3Df::faces_enum::top, BoundingBox3Df::faces_enum::bottom), + std::make_pair(BoundingBox3Df::faces_enum::front, BoundingBox3Df::faces_enum::back)}; // For each axis on which the feature touches both opposing periodic faces, the naive arithmetic centroid // lands in the empty middle of the wrapped feature. Replace that component with the minimum-image mean @@ -348,44 +348,4 @@ bool AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const Bo return isPeriodic; } - -bool AdjustCentroidsForPeriodicFaces(const ImageGeom& imageGeom, const UInt64AbstractDataStore& xRanges, const UInt64AbstractDataStore& yRanges, const UInt64AbstractDataStore& zRanges, - Float32AbstractDataStore& centroids) -{ - const size_t xPoints = imageGeom.getNumXCells() - 1; - const size_t yPoints = imageGeom.getNumYCells() - 1; - const size_t zPoints = imageGeom.getNumZCells() - 1; - - // Centroids are stored in physical coordinates (origin + (index + 0.5) * spacing), so the periodic - // wrap offset must also be a physical distance. Scale the (cell-count) span by the spacing on each - // axis; the offset is half the physical distance between the first and last cell centers. Using the - // raw cell count here would only be correct for unit spacing. - const auto spacing = imageGeom.getSpacing(); - - bool isAdjusted = false; - - const usize numFeatures = xRanges.size() / 2; - for(usize featureId = 0; featureId < numFeatures; featureId++) - { - if(xRanges[featureId * 2 + 0] == 0 && xRanges[featureId * 2 + 1] == xPoints) - { - isAdjusted = true; - centroids[featureId * 3 + 0] += (static_cast(xPoints) * spacing[0]) / 2.0f; - } - - if(yRanges[featureId * 2 + 0] == 0 && yRanges[featureId * 2 + 1] == yPoints) - { - isAdjusted = true; - centroids[featureId * 3 + 1] += (static_cast(yPoints) * spacing[1]) / 2.0f; - } - - if(zRanges[featureId * 2 + 0] == 0 && zRanges[featureId * 2 + 1] == zPoints) - { - isAdjusted = true; - centroids[featureId * 3 + 2] += (static_cast(zPoints) * spacing[2]) / 2.0f; - } - } - - return isAdjusted; -} } // namespace nx::core::GeometryHelpers::Topology diff --git a/src/simplnx/Utilities/GeometryHelpers.hpp b/src/simplnx/Utilities/GeometryHelpers.hpp index f8799e8b2f..0f5f565578 100644 --- a/src/simplnx/Utilities/GeometryHelpers.hpp +++ b/src/simplnx/Utilities/GeometryHelpers.hpp @@ -978,9 +978,6 @@ BoundingBoxFaces SIMPLNX_EXPORT FindElementPeriodicFaces(const BoundingBox3Df& b bool SIMPLNX_EXPORT AdjustCentroidsForPeriodicFaces(const BoundingBox3Df& boundingBox, const BoundingBoxFaces& faces, const Float32AbstractDataStore& vertices, const std::set& vertexSet, Float32AbstractDataStore& centroids, IGeometry::MeshIndexType featureId); -bool SIMPLNX_EXPORT AdjustCentroidsForPeriodicFaces(const ImageGeom& imageGeom, const UInt64AbstractDataStore& xRanges, const UInt64AbstractDataStore& yRanges, const UInt64AbstractDataStore& zRanges, - Float32AbstractDataStore& centroids); - /** * @brief * @tparam T