Merge pull request InsightSoftwareConsortium#6017 from hjmjohnson/fix-voronoi-infinite-loop

BUG: Fix infinite loop in VoronoiDiagram2DGenerator

Author: hjmjohnson

Modules/Segmentation/Voronoi/include/itkVoronoiDiagram2DGenerator.hxx

@@ -273,12 +273,56 @@ VoronoiDiagram2DGenerator<TCoordinate>::ConstructDiagram()
     buildEdges.push_back(curr);
     EdgeInfo front = curr;
     EdgeInfo back = curr;
-    while (!(rawEdges[i].empty()))
+    // Assemble raw edges into a connected chain for this Voronoi cell.
+    // Each iteration pops an edge from the deque and attempts to attach
+    // it to the front or back of the growing chain.  Edges that cannot
+    // attach are pushed back for retry, because later attachments change
+    // the chain endpoints and may make previously unattachable edges
+    // attachable.
+    //
+    // A stall counter tracks progress: it resets whenever an edge
+    // attaches, and terminates the loop when a full pass through the
+    // deque makes no progress.  Without this, certain degenerate seed
+    // configurations (near-collinear seeds, ITK issue #4386) cause an
+    // infinite loop because Fortune's algorithm produces near-zero-length
+    // edges whose endpoints differ by less than floating-point tolerance
+    // but have different vertex IDs.  These degenerate edges cannot
+    // attach because:
+    //   1. Their endpoints don't match any chain vertex by ID.
+    //   2. The chain may already be closed (front[0] == back[1]).
+    //   3. The boundary-bridging logic doesn't apply when the chain
+    //      endpoints are interior (not on the domain boundary).
+    // Such edges are safely dropped — they represent floating-point
+    // artifacts where two boundary intersection points should be
+    // identical in exact arithmetic.
+    auto remainingBeforeStall = rawEdges[i].size();
+    while (!(rawEdges[i].empty()) && (remainingBeforeStall != 0))
     {
+      --remainingBeforeStall;
       curr = rawEdges[i].front();
       rawEdges[i].pop_front();
+
+      // Check if this edge is a degenerate near-zero-length artifact.
+      // Fortune's algorithm can produce edges whose two endpoints map
+      // to the same geometric point (within DIFF_TOLERENCE) but have
+      // different vertex IDs.  These carry no geometric information
+      // and can be safely discarded.
+      const PointType & edgeStart = m_OutputVD->GetVertex(curr[0]);
+      const PointType & edgeEnd = m_OutputVD->GetVertex(curr[1]);
+      if (!differentPoint(edgeStart, edgeEnd))
+      {
+        itkDebugMacro("Dropping degenerate near-zero-length edge ["
+                      << curr[0] << " (" << edgeStart[0] << "," << edgeStart[1] << ") -> " << curr[1] << " ("
+                      << edgeEnd[0] << "," << edgeEnd[1] << ")]"
+                      << " for cell " << i << ": endpoints within DIFF_TOLERENCE=" << DIFF_TOLERENCE);
+        // Count as progress — this edge is resolved (discarded).
+        remainingBeforeStall = rawEdges[i].size();
+        continue;
+      }
+
       unsigned char frontbnd = Pointonbnd(front[0]);
       unsigned char backbnd = Pointonbnd(back[1]);
+      bool          edgeAttached = true;
       if (curr[0] == back[1])
       {
         buildEdges.push_back(curr);
@@ -353,12 +397,30 @@ VoronoiDiagram2DGenerator<TCoordinate>::ConstructDiagram()
         else
         {
           rawEdges[i].push_back(curr);
+          edgeAttached = false;
         }
       }
       else
       {
         rawEdges[i].push_back(curr);
+        edgeAttached = false;
       }
+      if (edgeAttached)
+      {
+        // Progress was made — chain endpoints changed, so previously
+        // unattachable edges may now be attachable.
+        remainingBeforeStall = rawEdges[i].size();
+      }
+    }
+    // After assembly, all edges for this cell should have been either
+    // attached to the chain or identified as degenerate artifacts.
+    // Any remaining edges indicate an unexpected algorithmic failure.
+    if (!rawEdges[i].empty())
+    {
+      itkExceptionMacro("VoronoiDiagram2DGenerator::ConstructDiagram: "
+                        << rawEdges[i].size() << " non-degenerate edge(s) could not be "
+                        << "assembled into cell " << i << " boundary chain. "
+                        << "This indicates an unexpected geometric configuration.");
     }
 
     curr = buildEdges.front();

Modules/Segmentation/Voronoi/test/CMakeLists.txt

@@ -4,11 +4,19 @@ set(
   itkVoronoiSegmentationImageFilterTest.cxx
   itkVoronoiSegmentationRGBImageFilterTest.cxx
   itkVoronoiDiagram2DTest.cxx
+  itkVoronoiDiagram2DInfiniteLoopTest.cxx
   itkVoronoiPartitioningImageFilterTest.cxx
 )
 
 createtestdriver(ITKVoronoi "${ITKVoronoi-Test_LIBRARIES}" "${ITKVoronoiTests}")
 
+itk_add_test(
+  NAME itkVoronoiDiagram2DInfiniteLoopTest
+  COMMAND
+    ITKVoronoiTestDriver
+    itkVoronoiDiagram2DInfiniteLoopTest
+)
+
 itk_add_test(
   NAME itkVoronoiSegmentationImageFilterTest
   COMMAND

Modules/Segmentation/Voronoi/test/itkVoronoiDiagram2DInfiniteLoopTest.cxx

@@ -0,0 +1,78 @@
+/*=========================================================================
+ *
+ *  Copyright NumFOCUS
+ *
+ *  Licensed under the Apache License, Version 2.0 (the "License");
+ *  you may not use this file except in compliance with the License.
+ *  You may obtain a copy of the License at
+ *
+ *         https://www.apache.org/licenses/LICENSE-2.0.txt
+ *
+ *  Unless required by applicable law or agreed to in writing, software
+ *  distributed under the License is distributed on an "AS IS" BASIS,
+ *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ *  See the License for the specific language governing permissions and
+ *  limitations under the License.
+ *
+ *=========================================================================*/
+
+#include "itkVoronoiDiagram2DGenerator.h"
+#include "itkTestingMacros.h"
+
+// Regression test for ITK issue #4386: near-collinear seed configurations
+// caused an infinite loop in VoronoiDiagram2DGenerator::ConstructDiagram().
+//
+// Root cause: Fortune's algorithm produces near-zero-length edges when
+// boundary intersection points coincide within floating-point tolerance.
+// These degenerate edges have different vertex IDs but geometrically
+// identical endpoints (within DIFF_TOLERENCE = 0.001).  They cannot be
+// attached to the growing boundary chain because:
+//   1. Their vertex IDs don't match any chain endpoint.
+//   2. The chain may already be closed (front == back).
+//   3. Boundary-bridging logic doesn't apply when chain endpoints are
+//      interior (not on the domain boundary).
+// The fix explicitly detects and drops these degenerate edges using the
+// existing differentPoint() tolerance check.
+int
+itkVoronoiDiagram2DInfiniteLoopTest(int argc, char * argv[])
+{
+  if (argc != 1)
+  {
+    std::cerr << "Usage: " << itkNameOfTestExecutableMacro(argv) << std::endl;
+    return EXIT_FAILURE;
+  }
+
+  using VoronoiDiagramType = itk::VoronoiDiagram2D<double>;
+  using VoronoiGeneratorType = itk::VoronoiDiagram2DGenerator<double>;
+  using PointType = VoronoiDiagramType::PointType;
+
+  // Six near-collinear seeds (x in [-1.40, -1.21]) that produce a
+  // degenerate Voronoi edge on the left domain boundary where two
+  // intersection points are ~0.00003 apart.
+  auto vg = VoronoiGeneratorType::New();
+  vg->SetOrigin(PointType{ { -1.61569, -1.76726 } });
+  vg->SetBoundary(PointType{ { 1.60174, 1.76345 } });
+  vg->AddOneSeed(PointType{ { -1.39649, 0.322212 } });
+  vg->AddOneSeed(PointType{ { -1.30128, 0.231786 } });
+  vg->AddOneSeed(PointType{ { -1.21509, 0.0515039 } });
+  vg->AddOneSeed(PointType{ { -1.22364, -0.030281 } });
+  vg->AddOneSeed(PointType{ { -1.22125, -0.120815 } });
+  vg->AddOneSeed(PointType{ { -1.25159, -0.23593 } });
+
+  // Without the fix, this call loops forever.  With the fix, the
+  // degenerate near-zero-length edge is detected and dropped, and
+  // the exception guard ensures no non-degenerate edges are lost.
+  ITK_TRY_EXPECT_NO_EXCEPTION(vg->Update());
+
+  // Verify all 6 cells were constructed with valid boundaries
+  auto vd = vg->GetOutput();
+  for (unsigned int i = 0; i < 6; ++i)
+  {
+    VoronoiDiagramType::CellAutoPointer cellPtr;
+    vd->GetCellId(i, cellPtr);
+    ITK_TEST_EXPECT_TRUE(cellPtr->GetNumberOfPoints() >= 2);
+  }
+
+  std::cout << "Test passed." << std::endl;
+  return EXIT_SUCCESS;
+}