From 2cb0c2fd675ee3091a4a029b0d321ced6d3d4176 Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:00 -0500
Subject: [PATCH 01/10] BUG: Guard N4 histogram sharpening against a degenerate
 bin range

A constant included region made histogramSlope zero, so the histogram
bin index became 0.0/0.0 = NaN whose float-to-int conversion is
undefined behavior feeding unchecked indexing. Independently, the
else-if min/max scan never assigned binMinimum when the traversed
values were strictly increasing.

Make the scans independent and degrade sharpening to an identity
mapping when the intensity range is degenerate.

Issue #6575, item B2.
---
 .../itkN4BiasFieldCorrectionImageFilter.hxx   | 12 ++++-
 .../BiasCorrection/test/CMakeLists.txt        |  6 ++-
 ...kN4BiasFieldCorrectionImageFilterGTest.cxx | 45 +++++++++++++++++++
 3 files changed, 61 insertions(+), 2 deletions(-)
 create mode 100644 Modules/Filtering/BiasCorrection/test/itkN4BiasFieldCorrectionImageFilterGTest.cxx

diff --git a/Modules/Filtering/BiasCorrection/include/itkN4BiasFieldCorrectionImageFilter.hxx b/Modules/Filtering/BiasCorrection/include/itkN4BiasFieldCorrectionImageFilter.hxx
index aa4fd7d6ae5..70729bb01fa 100644
--- a/Modules/Filtering/BiasCorrection/include/itkN4BiasFieldCorrectionImageFilter.hxx
+++ b/Modules/Filtering/BiasCorrection/include/itkN4BiasFieldCorrectionImageFilter.hxx
@@ -38,6 +38,7 @@ ITK_GCC_SUPPRESS_Wfloat_equal
 #include "vnl/vnl_complex_traits.h"
 #include "complex"
 #include "itkPrintHelper.h"
+#include <algorithm> // For clamp and copy.
 ITK_GCC_PRAGMA_POP
 
 namespace itk
@@ -285,12 +286,21 @@ N4BiasFieldCorrectionImageFilter<TInputImage, TMaskImage, TOutputImage>::Sharpen
       {
         binMaximum = pixel;
       }
-      else if (pixel < binMinimum)
+      if (pixel < binMinimum)
       {
         binMinimum = pixel;
       }
     }
   }
+  if (binMaximum <= binMinimum)
+  {
+    // Degenerate intensity range: sharpening is an identity mapping.
+    sharpenedImage->FillBuffer(0);
+    const ImageBufferRange identityImageBufferRange{ *sharpenedImage };
+    std::copy(
+      unsharpenedImageBufferRange.cbegin(), unsharpenedImageBufferRange.cend(), identityImageBufferRange.begin());
+    return;
+  }
   const RealType histogramSlope = (binMaximum - binMinimum) / static_cast<RealType>(this->m_NumberOfHistogramBins - 1);
 
   // Create the intensity profile (within the masked region, if applicable)
diff --git a/Modules/Filtering/BiasCorrection/test/CMakeLists.txt b/Modules/Filtering/BiasCorrection/test/CMakeLists.txt
index 96d7d3a8a45..7b525cece3e 100644
--- a/Modules/Filtering/BiasCorrection/test/CMakeLists.txt
+++ b/Modules/Filtering/BiasCorrection/test/CMakeLists.txt
@@ -7,7 +7,11 @@ set(
 
 createtestdriver(ITKBiasCorrection "${ITKBiasCorrection-Test_LIBRARIES}" "${ITKBiasCorrectionTests}")
 
-set(ITKBiasCorrectionGTests itkCompositeValleyFunctionGTest.cxx)
+set(
+  ITKBiasCorrectionGTests
+  itkCompositeValleyFunctionGTest.cxx
+  itkN4BiasFieldCorrectionImageFilterGTest.cxx
+)
 
 creategoogletestdriver(ITKBiasCorrection "${ITKBiasCorrection-Test_LIBRARIES}" "${ITKBiasCorrectionGTests}")
 
diff --git a/Modules/Filtering/BiasCorrection/test/itkN4BiasFieldCorrectionImageFilterGTest.cxx b/Modules/Filtering/BiasCorrection/test/itkN4BiasFieldCorrectionImageFilterGTest.cxx
new file mode 100644
index 00000000000..77df8e5f868
--- /dev/null
+++ b/Modules/Filtering/BiasCorrection/test/itkN4BiasFieldCorrectionImageFilterGTest.cxx
@@ -0,0 +1,45 @@
+/*=========================================================================
+ *
+ *  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 "itkN4BiasFieldCorrectionImageFilter.h"
+#include "itkImageRegionConstIterator.h"
+#include "itkGTest.h"
+
+// A constant image has a degenerate histogram range; sharpening must degrade to an
+// identity mapping instead of producing NaN bin indices (issue #6575, B2).
+TEST(N4BiasFieldCorrectionImageFilter, ConstantImageProducesFiniteOutput)
+{
+  using ImageType = itk::Image<float, 2>;
+  auto image = ImageType::New();
+  image->SetRegions(ImageType::RegionType{ ImageType::IndexType{}, ImageType::SizeType::Filled(32) });
+  image->Allocate();
+  image->FillBuffer(1.0F);
+
+  using FilterType = itk::N4BiasFieldCorrectionImageFilter<ImageType>;
+  auto filter = FilterType::New();
+  filter->SetInput(image);
+  filter->SetMaximumNumberOfIterations(FilterType::VariableSizeArrayType(1, 2));
+  filter->Update();
+
+  itk::ImageRegionConstIterator<ImageType> it(filter->GetOutput(), filter->GetOutput()->GetBufferedRegion());
+  for (it.GoToBegin(); !it.IsAtEnd(); ++it)
+  {
+    ASSERT_TRUE(std::isfinite(it.Get()));
+    ASSERT_NEAR(it.Get(), 1.0F, 1e-2) << " constant input must be (nearly) unchanged";
+  }
+}

From cd89b258283103d3d0030bbc32b9ab24783cc06a Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:04 -0500
Subject: [PATCH 02/10] BUG: Keep previous SLIC centroid for an empty cluster

A cluster ending an iteration with zero members divided its
accumulator by zero, leaving NaN spatial components whose later
float-to-integer rounding is undefined behavior feeding image reads.

Issue #6575, item B3.
---
 .../SuperPixel/include/itkSLICImageFilter.hxx       | 13 ++++++++++---
 1 file changed, 10 insertions(+), 3 deletions(-)

diff --git a/Modules/Segmentation/SuperPixel/include/itkSLICImageFilter.hxx b/Modules/Segmentation/SuperPixel/include/itkSLICImageFilter.hxx
index 443c02af2bd..576217ca9e4 100644
--- a/Modules/Segmentation/SuperPixel/include/itkSLICImageFilter.hxx
+++ b/Modules/Segmentation/SuperPixel/include/itkSLICImageFilter.hxx
@@ -642,10 +642,17 @@ SLICImageFilter<TInputImage, TOutputImage, TDistancePixel>::GenerateData()
     for (size_t i = 0; i * numberOfClusterComponents < m_Clusters.size(); ++i)
     {
 
-      ClusterType cluster(numberOfClusterComponents, &m_Clusters[i * numberOfClusterComponents]);
-      cluster /= clusterCount[i];
-
+      ClusterType       cluster(numberOfClusterComponents, &m_Clusters[i * numberOfClusterComponents]);
       const ClusterType oldCluster(numberOfClusterComponents, &m_OldClusters[i * numberOfClusterComponents]);
+      if (clusterCount[i] > 0)
+      {
+        cluster /= clusterCount[i];
+      }
+      else
+      {
+        // Empty cluster: keep the previous centroid to avoid a NaN center.
+        cluster.copy_in(oldCluster.data_block());
+      }
       l1Residual += Distance(cluster, oldCluster);
     }
 

From e0149edd35932242e234fbf9821362584e253241 Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:07 -0500
Subject: [PATCH 03/10] BUG: Guard Wiener deconvolution against zero
 signal-to-noise power

Pn / (Pf - Pn) divided by zero when the input power at a frequency
equaled the noise power spectral density constant; the magnitude
guard passes on the resulting infinity and NaN reached the output.
Return zero for such frequencies: zero estimated signal power means
nothing to recover.

Issue #6575, item B4.
---
 .../itkWienerDeconvolutionImageFilter.h       |  9 ++--
 .../Deconvolution/test/CMakeLists.txt         |  1 +
 ...itkWienerDeconvolutionImageFilterGTest.cxx | 41 +++++++++++++++++++
 3 files changed, 48 insertions(+), 3 deletions(-)
 create mode 100644 Modules/Filtering/Deconvolution/test/itkWienerDeconvolutionImageFilterGTest.cxx

diff --git a/Modules/Filtering/Deconvolution/include/itkWienerDeconvolutionImageFilter.h b/Modules/Filtering/Deconvolution/include/itkWienerDeconvolutionImageFilter.h
index 0c992790222..9494d46183b 100644
--- a/Modules/Filtering/Deconvolution/include/itkWienerDeconvolutionImageFilter.h
+++ b/Modules/Filtering/Deconvolution/include/itkWienerDeconvolutionImageFilter.h
@@ -170,11 +170,14 @@ class ITK_TEMPLATE_EXPORT WienerDeconvolutionFunctor
     // minus the power spectral density of the noise.
     TPixel Pf = std::norm(I);
 
-    TPixel denominator = std::norm(H) + (Pn / (Pf - Pn));
     TPixel value{};
-    if (itk::Math::Absolute(denominator) >= m_KernelZeroMagnitudeThreshold)
+    if (Pf != Pn)
     {
-      value = I * (std::conj(H) / denominator);
+      const TPixel denominator = std::norm(H) + (Pn / (Pf - Pn));
+      if (itk::Math::Absolute(denominator) >= m_KernelZeroMagnitudeThreshold)
+      {
+        value = I * (std::conj(H) / denominator);
+      }
     }
 
     return value;
diff --git a/Modules/Filtering/Deconvolution/test/CMakeLists.txt b/Modules/Filtering/Deconvolution/test/CMakeLists.txt
index bb8e3ed6319..25aa381788a 100644
--- a/Modules/Filtering/Deconvolution/test/CMakeLists.txt
+++ b/Modules/Filtering/Deconvolution/test/CMakeLists.txt
@@ -15,6 +15,7 @@ createtestdriver(ITKDeconvolution "${ITKDeconvolution-Test_LIBRARIES}" "${ITKDec
 set(
   ITKDeconvolutionGTests
   itkProjectedIterativeDeconvolutionImageFilterGTest.cxx
+  itkWienerDeconvolutionImageFilterGTest.cxx
 )
 
 creategoogletestdriver(ITKDeconvolution "${ITKDeconvolution-Test_LIBRARIES}" "${ITKDeconvolutionGTests}")
diff --git a/Modules/Filtering/Deconvolution/test/itkWienerDeconvolutionImageFilterGTest.cxx b/Modules/Filtering/Deconvolution/test/itkWienerDeconvolutionImageFilterGTest.cxx
new file mode 100644
index 00000000000..0337734cc6c
--- /dev/null
+++ b/Modules/Filtering/Deconvolution/test/itkWienerDeconvolutionImageFilterGTest.cxx
@@ -0,0 +1,41 @@
+/*=========================================================================
+ *
+ *  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 "itkWienerDeconvolutionImageFilter.h"
+#include "itkGTest.h"
+#include <complex>
+
+// When the input power at a frequency equals the noise power spectral density, the
+// regularization term divides by zero; the functor must return zero instead of NaN
+// (issue #6575, B4).
+TEST(WienerDeconvolutionFunctor, ZeroSignalToNoiseYieldsZero)
+{
+  using ComplexType = std::complex<double>;
+  itk::Functor::WienerDeconvolutionFunctor<ComplexType> functor;
+  functor.SetNoisePowerSpectralDensityConstant(4.0);
+  functor.SetKernelZeroMagnitudeThreshold(1.0e-4);
+
+  // |I|^2 == Pn triggers the degenerate branch.
+  const ComplexType value = functor(ComplexType(2.0, 0.0), ComplexType(1.0, 0.0));
+  EXPECT_EQ(value, ComplexType(0.0, 0.0));
+
+  // A non-degenerate frequency still deconvolves finitely.
+  const ComplexType regular = functor(ComplexType(3.0, 0.0), ComplexType(1.0, 0.0));
+  EXPECT_TRUE(std::isfinite(regular.real()) && std::isfinite(regular.imag()));
+  EXPECT_NE(regular, ComplexType(0.0, 0.0));
+}

From 0c70d47d3c6390366a0b86328922ac85fd4e3476 Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:30 -0500
Subject: [PATCH 04/10] BUG: Make adaptive histogram equalization of a constant
 image an identity

A constant image gave a zero intensity range, so the normalization
divided by zero and every output pixel became NaN.

Issue #6575, item B9.
---
 .../itkAdaptiveEqualizationHistogram.h        |  5 +++
 .../ImageStatistics/test/CMakeLists.txt       |  1 +
 ...eHistogramEqualizationImageFilterGTest.cxx | 44 +++++++++++++++++++
 3 files changed, 50 insertions(+)
 create mode 100644 Modules/Filtering/ImageStatistics/test/itkAdaptiveHistogramEqualizationImageFilterGTest.cxx

diff --git a/Modules/Filtering/ImageStatistics/include/itkAdaptiveEqualizationHistogram.h b/Modules/Filtering/ImageStatistics/include/itkAdaptiveEqualizationHistogram.h
index fd84726eb19..7c587a1a6a4 100644
--- a/Modules/Filtering/ImageStatistics/include/itkAdaptiveEqualizationHistogram.h
+++ b/Modules/Filtering/ImageStatistics/include/itkAdaptiveEqualizationHistogram.h
@@ -74,6 +74,11 @@ class AdaptiveEqualizationHistogram
     // AdaptiveHistogramEqualization compute kernel components with
     // float, but use double for accumulate and temporaries.
     const double iscale = static_cast<double>(m_Maximum) - m_Minimum;
+    if (iscale == 0.0)
+    {
+      // Constant image: equalization is an identity mapping.
+      return static_cast<TOutputPixel>(pixel);
+    }
 
     double         sum = 0.0;
     auto           itMap = m_Map.begin();
diff --git a/Modules/Filtering/ImageStatistics/test/CMakeLists.txt b/Modules/Filtering/ImageStatistics/test/CMakeLists.txt
index 41c49a34924..5eb1da61d20 100644
--- a/Modules/Filtering/ImageStatistics/test/CMakeLists.txt
+++ b/Modules/Filtering/ImageStatistics/test/CMakeLists.txt
@@ -467,6 +467,7 @@ itk_add_test(
 
 set(
   ITKImageStatisticsGTests
+  itkAdaptiveHistogramEqualizationImageFilterGTest.cxx
   itkLabelOverlapMeasuresImageFilterGTest.cxx
   itkMinimumMaximumImageFilterGTest.cxx
 )
diff --git a/Modules/Filtering/ImageStatistics/test/itkAdaptiveHistogramEqualizationImageFilterGTest.cxx b/Modules/Filtering/ImageStatistics/test/itkAdaptiveHistogramEqualizationImageFilterGTest.cxx
new file mode 100644
index 00000000000..20b55421576
--- /dev/null
+++ b/Modules/Filtering/ImageStatistics/test/itkAdaptiveHistogramEqualizationImageFilterGTest.cxx
@@ -0,0 +1,44 @@
+/*=========================================================================
+ *
+ *  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 "itkAdaptiveHistogramEqualizationImageFilter.h"
+#include "itkImageRegionConstIterator.h"
+#include "itkGTest.h"
+
+// A constant image has zero intensity range; equalization must be an identity mapping
+// instead of dividing by zero and producing NaN everywhere (issue #6575, B9).
+TEST(AdaptiveHistogramEqualizationImageFilter, ConstantImageIsIdentity)
+{
+  using ImageType = itk::Image<float, 2>;
+  auto image = ImageType::New();
+  image->SetRegions(ImageType::RegionType{ ImageType::IndexType{}, ImageType::SizeType::Filled(16) });
+  image->Allocate();
+  image->FillBuffer(7.0F);
+
+  using FilterType = itk::AdaptiveHistogramEqualizationImageFilter<ImageType>;
+  auto filter = FilterType::New();
+  filter->SetInput(image);
+  filter->SetRadius(2);
+  filter->Update();
+
+  itk::ImageRegionConstIterator<ImageType> it(filter->GetOutput(), filter->GetOutput()->GetBufferedRegion());
+  for (it.GoToBegin(); !it.IsAtEnd(); ++it)
+  {
+    ASSERT_EQ(it.Get(), 7.0F);
+  }
+}

From cde512cc4c45234a46fcf103ac22e5139d2317b7 Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:36 -0500
Subject: [PATCH 05/10] BUG: Guard STAPLE estimates against a zero weight
 denominator

All-background inputs zero the weight image, making the sensitivity
denominator 0.0/0.0 = NaN with no guard. Initialize the estimates and
keep the previous value when a denominator vanishes.

Issue #6575, item B11.
---
 .../include/itkSTAPLEImageFilter.hxx          | 12 +++-
 .../ImageCompare/test/CMakeLists.txt          |  3 +
 .../test/itkSTAPLEImageFilterGTest.cxx        | 56 +++++++++++++++++++
 3 files changed, 69 insertions(+), 2 deletions(-)
 create mode 100644 Modules/Filtering/ImageCompare/test/itkSTAPLEImageFilterGTest.cxx

diff --git a/Modules/Filtering/ImageCompare/include/itkSTAPLEImageFilter.hxx b/Modules/Filtering/ImageCompare/include/itkSTAPLEImageFilter.hxx
index 844df742b33..90647ea6a65 100644
--- a/Modules/Filtering/ImageCompare/include/itkSTAPLEImageFilter.hxx
+++ b/Modules/Filtering/ImageCompare/include/itkSTAPLEImageFilter.hxx
@@ -73,6 +73,8 @@ STAPLEImageFilter<TInputImage, TOutputImage>::GenerateData()
   {
     last_p[i] = -10.0;
     last_q[i] = -10.0;
+    p[i] = 1.0;
+    q[i] = 1.0;
   }
 
   // Come up with an initial Wi which is simply the average of
@@ -158,8 +160,14 @@ STAPLEImageFilter<TInputImage, TOutputImage>::GenerateData()
         out.NextLine();
       }
 
-      p[i] = p_num / p_denom;
-      q[i] = q_num / q_denom;
+      if (p_denom != 0.0)
+      {
+        p[i] = p_num / p_denom;
+      }
+      if (q_denom != 0.0)
+      {
+        q[i] = q_num / q_denom;
+      }
     }
 
     // Now recreate W using the new p's and q's
diff --git a/Modules/Filtering/ImageCompare/test/CMakeLists.txt b/Modules/Filtering/ImageCompare/test/CMakeLists.txt
index 87260430bd1..a7dc12111f3 100644
--- a/Modules/Filtering/ImageCompare/test/CMakeLists.txt
+++ b/Modules/Filtering/ImageCompare/test/CMakeLists.txt
@@ -13,6 +13,9 @@ set(
 
 createtestdriver(ITKImageCompare "${ITKImageCompare-Test_LIBRARIES}" "${ITKImageCompareTests}")
 
+set(ITKImageCompareGTests itkSTAPLEImageFilterGTest.cxx)
+creategoogletestdriver(ITKImageCompare "${ITKImageCompare-Test_LIBRARIES}" "${ITKImageCompareGTests}")
+
 itk_add_test(
   NAME itkAbsoluteValueDifferenceImageFilterTest
   COMMAND
diff --git a/Modules/Filtering/ImageCompare/test/itkSTAPLEImageFilterGTest.cxx b/Modules/Filtering/ImageCompare/test/itkSTAPLEImageFilterGTest.cxx
new file mode 100644
index 00000000000..5d4e4fd7c58
--- /dev/null
+++ b/Modules/Filtering/ImageCompare/test/itkSTAPLEImageFilterGTest.cxx
@@ -0,0 +1,56 @@
+/*=========================================================================
+ *
+ *  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 "itkSTAPLEImageFilter.h"
+#include "itkImageRegionConstIterator.h"
+#include "itkGTest.h"
+
+// All-background segmentations zero the weight image, so the sensitivity denominator is
+// zero; the estimates must stay finite instead of becoming NaN (issue #6575, B11).
+TEST(STAPLEImageFilter, AllBackgroundInputsStayFinite)
+{
+  using InputImageType = itk::Image<unsigned char, 2>;
+  using OutputImageType = itk::Image<double, 2>;
+
+  const auto makeImage = [] {
+    auto image = InputImageType::New();
+    image->SetRegions(InputImageType::RegionType{ InputImageType::IndexType{}, InputImageType::SizeType::Filled(16) });
+    image->Allocate();
+    image->FillBuffer(0);
+    return image;
+  };
+
+  using FilterType = itk::STAPLEImageFilter<InputImageType, OutputImageType>;
+  auto filter = FilterType::New();
+  filter->SetInput(0, makeImage());
+  filter->SetInput(1, makeImage());
+  filter->SetForegroundValue(1);
+  filter->SetMaximumIterations(5);
+  filter->Update();
+
+  for (unsigned int i = 0; i < 2; ++i)
+  {
+    EXPECT_TRUE(std::isfinite(filter->GetSensitivity(i))) << " input " << i;
+    EXPECT_TRUE(std::isfinite(filter->GetSpecificity(i))) << " input " << i;
+  }
+  itk::ImageRegionConstIterator<OutputImageType> it(filter->GetOutput(), filter->GetOutput()->GetBufferedRegion());
+  for (it.GoToBegin(); !it.IsAtEnd(); ++it)
+  {
+    ASSERT_TRUE(std::isfinite(it.Get()));
+  }
+}

From 50741405d253a1a41aa2fffabacecf4d701c084e Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:40 -0500
Subject: [PATCH 06/10] BUG: Test the false-positive-error denominator, not the
 source count

The degenerate-value guard tested the source voxel count while the
division is by FP + TN, so a label covering the whole image passed the
guard and produced 0.0/0.0 = NaN.

Issue #6575, item B12.
---
 .../itkLabelOverlapMeasuresImageFilter.hxx       | 11 +++++------
 .../itkLabelOverlapMeasuresImageFilterGTest.cxx  | 16 ++++++++++++++++
 2 files changed, 21 insertions(+), 6 deletions(-)

diff --git a/Modules/Filtering/ImageStatistics/include/itkLabelOverlapMeasuresImageFilter.hxx b/Modules/Filtering/ImageStatistics/include/itkLabelOverlapMeasuresImageFilter.hxx
index bbcac5cbdb1..a565f4d2dc7 100644
--- a/Modules/Filtering/ImageStatistics/include/itkLabelOverlapMeasuresImageFilter.hxx
+++ b/Modules/Filtering/ImageStatistics/include/itkLabelOverlapMeasuresImageFilter.hxx
@@ -385,17 +385,16 @@ LabelOverlapMeasuresImageFilter<TLabelImage>::GetFalsePositiveError(LabelType la
     return 0.0;
   }
 
-  RealType value;
-  if (Math::ExactlyEquals(mapIt->second.m_Source, 0.0))
+  RealType   value;
+  const auto nComplementIntersection = nVox - mapIt->second.m_Union; // TN
+  const auto denominator = mapIt->second.m_SourceComplement + nComplementIntersection;
+  if (denominator == 0)
   {
     value = NumericTraits<RealType>::max();
   }
   else
   {
-    auto nComplementIntersection = nVox - mapIt->second.m_Union; // TN
-
-    value = static_cast<RealType>(mapIt->second.m_SourceComplement) /
-            static_cast<RealType>(mapIt->second.m_SourceComplement + nComplementIntersection);
+    value = static_cast<RealType>(mapIt->second.m_SourceComplement) / static_cast<RealType>(denominator);
   }
 
   return value;
diff --git a/Modules/Filtering/ImageStatistics/test/itkLabelOverlapMeasuresImageFilterGTest.cxx b/Modules/Filtering/ImageStatistics/test/itkLabelOverlapMeasuresImageFilterGTest.cxx
index eddd1960e08..e2638d7f571 100644
--- a/Modules/Filtering/ImageStatistics/test/itkLabelOverlapMeasuresImageFilterGTest.cxx
+++ b/Modules/Filtering/ImageStatistics/test/itkLabelOverlapMeasuresImageFilterGTest.cxx
@@ -204,3 +204,19 @@ TEST_F(LabelOverlapMeasuresImageFilterFixture, test3)
   EXPECT_NEAR(filter->GetFalsePositiveError(2), 4.76837612950e-07, 1e-17);
   EXPECT_NEAR(filter->GetFalseDiscoveryRate(2), 1.0 / 3.0, 0.0);
 }
+
+// A label covering the whole image has FP = TN = 0; the zero denominator must map to
+// the degenerate-value convention, not NaN (issue #6575, B12).
+TEST_F(LabelOverlapMeasuresImageFilterFixture, WholeImageLabelFalsePositiveError)
+{
+  auto source = Utils::CreateImage(1);
+  auto target = Utils::CreateImage(1);
+
+  auto filter = Utils::FilterType::New();
+  filter->SetSourceImage(source);
+  filter->SetTargetImage(target);
+  filter->Update();
+
+  using RealType = Utils::FilterType::RealType;
+  EXPECT_EQ(filter->GetFalsePositiveError(1), itk::NumericTraits<RealType>::max());
+}

From 75e4fdad838ed1b54af23c33f96f5579f03467f0 Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:45 -0500
Subject: [PATCH 07/10] BUG: Return zero sigma for a single-sample box window

A radius-0 window has one sample, so the unbiased variance divided by
zero and every output pixel became sqrt(0.0/0.0) = NaN.

Issue #6575, item B20.
---
 .../Smoothing/include/itkBoxUtilities.h       |  9 ++--
 .../Filtering/Smoothing/test/CMakeLists.txt   |  1 +
 .../test/itkBoxSigmaImageFilterGTest.cxx      | 47 +++++++++++++++++++
 3 files changed, 54 insertions(+), 3 deletions(-)
 create mode 100644 Modules/Filtering/Smoothing/test/itkBoxSigmaImageFilterGTest.cxx

diff --git a/Modules/Filtering/Smoothing/include/itkBoxUtilities.h b/Modules/Filtering/Smoothing/include/itkBoxUtilities.h
index e7476b0679e..151c49d3a54 100644
--- a/Modules/Filtering/Smoothing/include/itkBoxUtilities.h
+++ b/Modules/Filtering/Smoothing/include/itkBoxUtilities.h
@@ -473,7 +473,9 @@ BoxSigmaCalculatorFunction(const TInputImage *               accImage,
           ++(cornerItVec[k]);
         }
 
-        oIt.Set(static_cast<OutputPixelType>(std::sqrt((squareSum - sum * sum / pixelscount) / (pixelscount - 1))));
+        oIt.Set(pixelscount > 1
+                  ? static_cast<OutputPixelType>(std::sqrt((squareSum - sum * sum / pixelscount) / (pixelscount - 1)))
+                  : OutputPixelType{});
       }
     }
     else
@@ -538,8 +540,9 @@ BoxSigmaCalculatorFunction(const TInputImage *               accImage,
           }
         }
 
-        oIt.Set(
-          static_cast<OutputPixelType>(std::sqrt((squareSum - sum * sum / edgepixelscount) / (edgepixelscount - 1))));
+        oIt.Set(edgepixelscount > 1 ? static_cast<OutputPixelType>(
+                                        std::sqrt((squareSum - sum * sum / edgepixelscount) / (edgepixelscount - 1)))
+                                    : OutputPixelType{});
       }
     }
   }
diff --git a/Modules/Filtering/Smoothing/test/CMakeLists.txt b/Modules/Filtering/Smoothing/test/CMakeLists.txt
index d3864da347e..19f3cd69ba6 100644
--- a/Modules/Filtering/Smoothing/test/CMakeLists.txt
+++ b/Modules/Filtering/Smoothing/test/CMakeLists.txt
@@ -396,6 +396,7 @@ itk_add_test(
 
 set(
   ITKSmoothingGTests
+  itkBoxSigmaImageFilterGTest.cxx
   itkMeanImageFilterGTest.cxx
   itkMedianImageFilterGTest.cxx
 )
diff --git a/Modules/Filtering/Smoothing/test/itkBoxSigmaImageFilterGTest.cxx b/Modules/Filtering/Smoothing/test/itkBoxSigmaImageFilterGTest.cxx
new file mode 100644
index 00000000000..7746e833d55
--- /dev/null
+++ b/Modules/Filtering/Smoothing/test/itkBoxSigmaImageFilterGTest.cxx
@@ -0,0 +1,47 @@
+/*=========================================================================
+ *
+ *  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 "itkBoxSigmaImageFilter.h"
+#include "itkImageRegionConstIteratorWithIndex.h"
+#include "itkGTest.h"
+
+// A radius-0 window holds a single sample, so the (unbiased) standard deviation must be
+// zero, not sqrt(0/0) = NaN (issue #6575, B20).
+TEST(BoxSigmaImageFilter, RadiusZeroYieldsZeroNotNaN)
+{
+  using ImageType = itk::Image<float, 2>;
+  auto image = ImageType::New();
+  image->SetRegions(ImageType::RegionType{ ImageType::IndexType{}, ImageType::SizeType::Filled(8) });
+  image->Allocate();
+  for (itk::ImageRegionIteratorWithIndex<ImageType> it(image, image->GetBufferedRegion()); !it.IsAtEnd(); ++it)
+  {
+    it.Set(static_cast<float>(it.GetIndex()[0] + 10 * it.GetIndex()[1]));
+  }
+
+  using FilterType = itk::BoxSigmaImageFilter<ImageType, ImageType>;
+  auto filter = FilterType::New();
+  filter->SetInput(image);
+  filter->SetRadius(0);
+  filter->Update();
+
+  itk::ImageRegionConstIteratorWithIndex<ImageType> it(filter->GetOutput(), filter->GetOutput()->GetBufferedRegion());
+  for (it.GoToBegin(); !it.IsAtEnd(); ++it)
+  {
+    EXPECT_EQ(it.Get(), 0.0F) << " at index " << it.GetIndex();
+  }
+}

From bf581031dc17c6bc475636174e702bba3bdd154c Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:49 -0500
Subject: [PATCH 08/10] BUG: Guard structure-tensor rescaling against a zero
 maximum trace

A constant input has zero tensors, so the unit-maximum-trace rescale
divided by zero and multiplied the tensors into NaN, later driving an
undefined float-to-int conversion in the diffusion step count.

Issue #6575, items B23 and part of the B23 chain.
---
 .../include/itkStructureTensorImageFilter.hxx |  3 +-
 .../test/CMakeLists.txt                       |  1 +
 .../itkStructureTensorImageFilterGTest.cxx    | 51 +++++++++++++++++++
 3 files changed, 54 insertions(+), 1 deletion(-)
 create mode 100644 Modules/Filtering/AnisotropicDiffusionLBR/test/itkStructureTensorImageFilterGTest.cxx

diff --git a/Modules/Filtering/AnisotropicDiffusionLBR/include/itkStructureTensorImageFilter.hxx b/Modules/Filtering/AnisotropicDiffusionLBR/include/itkStructureTensorImageFilter.hxx
index b6d1cfaa77d..bf3e6389668 100644
--- a/Modules/Filtering/AnisotropicDiffusionLBR/include/itkStructureTensorImageFilter.hxx
+++ b/Modules/Filtering/AnisotropicDiffusionLBR/include/itkStructureTensorImageFilter.hxx
@@ -150,7 +150,8 @@ StructureTensorImageFilter<TImage, TTensorImage>::GenerateData()
 
   traceFilter->Update();
   maximumCalculator->ComputeMaximum();
-  m_PostRescaling = 1. / maximumCalculator->GetMaximum();
+  const auto maximum = maximumCalculator->GetMaximum();
+  m_PostRescaling = (maximum > 0.) ? 1. / maximum : 1.;
   scaleFilter->GetFunctor().scaling = m_PostRescaling;
   scaleFilter->Update();
   this->GraftOutput(scaleFilter->GetOutput());
diff --git a/Modules/Filtering/AnisotropicDiffusionLBR/test/CMakeLists.txt b/Modules/Filtering/AnisotropicDiffusionLBR/test/CMakeLists.txt
index 6be49446c8d..c47d032e291 100644
--- a/Modules/Filtering/AnisotropicDiffusionLBR/test/CMakeLists.txt
+++ b/Modules/Filtering/AnisotropicDiffusionLBR/test/CMakeLists.txt
@@ -7,6 +7,7 @@ createtestdriver(AnisotropicDiffusionLBR  "${AnisotropicDiffusionLBR-Test_LIBRAR
 set(
   AnisotropicDiffusionLBRGTests
   itkLinearAnisotropicDiffusionLBRImageFilterGTest.cxx
+  itkStructureTensorImageFilterGTest.cxx
 )
 
 creategoogletestdriver(AnisotropicDiffusionLBR "${AnisotropicDiffusionLBR-Test_LIBRARIES}" "${AnisotropicDiffusionLBRGTests}")
diff --git a/Modules/Filtering/AnisotropicDiffusionLBR/test/itkStructureTensorImageFilterGTest.cxx b/Modules/Filtering/AnisotropicDiffusionLBR/test/itkStructureTensorImageFilterGTest.cxx
new file mode 100644
index 00000000000..dac2ac7f429
--- /dev/null
+++ b/Modules/Filtering/AnisotropicDiffusionLBR/test/itkStructureTensorImageFilterGTest.cxx
@@ -0,0 +1,51 @@
+/*=========================================================================
+ *
+ *  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 "itkStructureTensorImageFilter.h"
+#include "itkImageRegionConstIterator.h"
+#include "itkGTest.h"
+
+// A constant image has zero gradients, so the tensor-trace maximum is zero; the
+// adimensionizing rescale must not divide by it and poison the tensors with NaN
+// (issue #6575, B23).
+TEST(StructureTensorImageFilter, ConstantImageYieldsFiniteZeroTensors)
+{
+  using ImageType = itk::Image<float, 2>;
+  auto image = ImageType::New();
+  image->SetRegions(ImageType::RegionType{ ImageType::IndexType{}, ImageType::SizeType::Filled(16) });
+  image->Allocate();
+  image->FillBuffer(0.0F);
+
+  using FilterType = itk::StructureTensorImageFilter<ImageType>;
+  auto filter = FilterType::New();
+  filter->SetInput(image);
+  filter->SetRescaleForUnitMaximumTrace(true);
+  filter->Update();
+
+  using TensorImageType = FilterType::TensorImageType;
+  itk::ImageRegionConstIterator<TensorImageType> it(filter->GetOutput(), filter->GetOutput()->GetBufferedRegion());
+  for (it.GoToBegin(); !it.IsAtEnd(); ++it)
+  {
+    const auto & tensor = it.Get();
+    for (unsigned int c = 0; c < tensor.Size(); ++c)
+    {
+      ASSERT_TRUE(std::isfinite(tensor[c]));
+      ASSERT_NEAR(tensor[c], 0.0F, 1e-20);
+    }
+  }
+}

From 86bb09e13791f9d6392a019b94a09197a152967f Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:53 -0500
Subject: [PATCH 09/10] BUG: Clamp denoising update to the representable pixel
 range

The GAUSSIAN fidelity branch and the pure-smoothing path can drive the
update negative; casting that to an unsigned pixel type is undefined
behavior.

Issue #6575, item B27.
---
 .../include/itkPatchBasedDenoisingImageFilter.hxx        | 9 +++++++++
 1 file changed, 9 insertions(+)

diff --git a/Modules/Filtering/Denoising/include/itkPatchBasedDenoisingImageFilter.hxx b/Modules/Filtering/Denoising/include/itkPatchBasedDenoisingImageFilter.hxx
index a929aa42d16..bb2b1ef4b67 100644
--- a/Modules/Filtering/Denoising/include/itkPatchBasedDenoisingImageFilter.hxx
+++ b/Modules/Filtering/Denoising/include/itkPatchBasedDenoisingImageFilter.hxx
@@ -30,6 +30,7 @@
 #include "itkSpatialNeighborSubsampler.h"
 #include "itkMacro.h"
 #include "itkMath.h"
+#include <algorithm> // For clamp and copy.
 
 namespace itk
 {
@@ -2075,6 +2076,14 @@ PatchBasedDenoisingImageFilter<TInputImage, TOutputImage>::ThreadedComputeImageU
       // this->ApplyUpdate() will then copy the results to outputIt since we
       // have to wait until all threads have finished using outputIt to
       // compute the update
+      for (unsigned int pc = 0; pc < m_NumPixelComponents; ++pc)
+      {
+        // Clamp to the representable range: the float-to-integer conversion of an
+        // out-of-range update is undefined behavior.
+        const auto lowest = static_cast<RealValueType>(NumericTraits<PixelValueType>::NonpositiveMin());
+        const auto highest = static_cast<RealValueType>(NumericTraits<PixelValueType>::max());
+        this->SetComponent(result, pc, std::clamp(this->GetComponent(result, pc), lowest, highest));
+      }
       updateIt.Set(static_cast<PixelType>(result));
 
       ++updateIt;

From 97e6a1169f0bcf5d05fc9250181f5dfc5efd3bb5 Mon Sep 17 00:00:00 2001
From: "Hans J. Johnson" <hans-johnson@uiowa.edu>
Date: Fri, 10 Jul 2026 20:47:56 -0500
Subject: [PATCH 10/10] BUG: Reject out-of-bound Gaussian variates before
 integer conversion

A normal variate is unbounded below and the float-to-unsigned
conversion of a negative value is undefined behavior; the rejection
loop relied on the wrapped result.

Issue #6575, item B29.
---
 .../itkGaussianRandomSpatialNeighborSubsampler.hxx   | 12 ++++++------
 1 file changed, 6 insertions(+), 6 deletions(-)

diff --git a/Modules/Numerics/Statistics/include/itkGaussianRandomSpatialNeighborSubsampler.hxx b/Modules/Numerics/Statistics/include/itkGaussianRandomSpatialNeighborSubsampler.hxx
index 987d519963b..44a47daeb8d 100644
--- a/Modules/Numerics/Statistics/include/itkGaussianRandomSpatialNeighborSubsampler.hxx
+++ b/Modules/Numerics/Statistics/include/itkGaussianRandomSpatialNeighborSubsampler.hxx
@@ -48,14 +48,14 @@ GaussianRandomSpatialNeighborSubsampler<TSample, TRegion>::GetIntegerVariate(Ran
     itkExceptionMacro("upperBound (" << upperBound << ") not >= to lowerBound(" << lowerBound << ')');
   }
 
-  RandomIntType randInt = 0;
-
+  // Reject out-of-bound variates before the integer conversion: converting a
+  // negative value to an unsigned type is undefined behavior.
+  RealType floorVar = 0;
   do
   {
-    const RealType randVar = this->m_RandomNumberGenerator->GetNormalVariate(mean, m_Variance);
-    randInt = static_cast<RandomIntType>(std::floor(randVar));
-  } while ((randInt < lowerBound) || (randInt > upperBound));
-  return randInt;
+    floorVar = std::floor(this->m_RandomNumberGenerator->GetNormalVariate(mean, m_Variance));
+  } while ((floorVar < static_cast<RealType>(lowerBound)) || (floorVar > static_cast<RealType>(upperBound)));
+  return static_cast<RandomIntType>(floorVar);
 }
 
 template <typename TSample, typename TRegion>