Repository: InsightSoftwareConsortium/ITK (language: python)

Issue:
`PointSet::SetPoints(PointsVectorContainer *)` overload leads to undefined behavior
### Description

`PointSet::SetPoints(PointsVectorContainer * points)`  internally casts its argument to a `PointsContainer *`:

https://github.com/InsightSoftwareConsortium/ITK/blob/9a5dade7f7e13fcef307f8224a01c2ce2b0ab13b/Modules/Core/Common/include/itkPointSet.hxx#L76

This leads to undefined behavior, when the `PointsVectorContainer` object is being used as a `PointsContainer` object.

### Steps to Reproduce

The test code already has undefined behavior, even if it might just work on the currently tested platforms:

https://github.com/InsightSoftwareConsortium/ITK/blob/9a5dade7f7e13fcef307f8224a01c2ce2b0ab13b/Modules/Core/Common/test/itkPointSetTest.cxx#L102-L127

### Versions

- The `PointSet::SetPoints(PointsVectorContainer *)` overload  is introduced by pull request #3154 commit 6827360fa6a06fbe073911d7f8554c2c3368721d "ENH: Adding SetPoints in PointSet for 1D input array", merged on Mar 16, 2022, and included with both [v5.3.0](https://github.com/InsightSoftwareConsortium/ITK/releases/tag/v5.3.0) and [v5.4.0](https://github.com/InsightSoftwareConsortium/ITK/releases/tag/v5.4.0).


### Additional Information

Discussed before at https://github.com/InsightSoftwareConsortium/ITK/pull/3154#discussion_r1755213148, with @PranjalSahu


Discussion / hints:

Fix patch:
diff --git a/Modules/Core/Common/include/itkPointSet.h b/Modules/Core/Common/include/itkPointSet.h
index 9b8f6d75e94..9300bc72488 100644
--- a/Modules/Core/Common/include/itkPointSet.h
+++ b/Modules/Core/Common/include/itkPointSet.h
@@ -164,10 +164,15 @@ class ITK_TEMPLATE_EXPORT PointSet : public DataObject
   void
   SetPoints(PointsContainer *);
 
-  /** Set the points container using a 1D vector. */
+  /** Set the points container using a 1D vector.
+  \warning This member function is unsafe. It may just work, but it may also lead to undefined behavior. */
   void
   SetPoints(PointsVectorContainer *);
 
+  /** Sets the points by specifying its coordinates. */
+  void
+  SetPointsByCoordinates(const std::vector<CoordRepType> & coordinates);
+
   /** Get the points container. */
   PointsContainer *
   GetPoints();
diff --git a/Modules/Core/Common/include/itkPointSet.hxx b/Modules/Core/Common/include/itkPointSet.hxx
index 0aed6b93cae..c868dfc01a6 100644
--- a/Modules/Core/Common/include/itkPointSet.hxx
+++ b/Modules/Core/Common/include/itkPointSet.hxx
@@ -73,12 +73,64 @@ PointSet<TPixelType, VDimension, TMeshTraits>::SetPoints(PointsVectorContainer *
   {
     itkExceptionMacro("Number of entries in given 1d array incompatible with the point dimension");
   }
+
+  // Note: this cast is unsafe. It may lead to undefined behavior.
   auto * pointsPtr = reinterpret_cast<PointsContainer *>(points);
 
   m_PointsContainer = pointsPtr;
   this->Modified();
 }
 
+
+template <typename TPixelType, unsigned int VDimension, typename TMeshTraits>
+void
+PointSet<TPixelType, VDimension, TMeshTraits>::SetPointsByCoordinates(const std::vector<CoordRepType> & coordinates)
+{
+  itkDebugMacro("Setting the points to the specified coordinates");
+
+  const size_t numberOfCoordinates = coordinates.size();
+
+  if (numberOfCoordinates % PointDimension != 0)
+  {
+    itkExceptionMacro("Number of specified coordinates incompatible with the point dimension");
+  }
+
+  const size_t numberOfPoints = numberOfCoordinates / PointDimension;
+
+  if (m_PointsContainer == nullptr)
+  {
+    m_PointsContainer = PointsContainer::New();
+  }
+
+  using STLContainerType = typename PointsContainer::STLContainerType;
+
+  STLContainerType & points = m_PointsContainer->CastToSTLContainer();
+  points.clear();
+
+  if constexpr (std::is_same_v<STLContainerType, std::vector<PointType>>)
+  {
+    // STLContainerType is either an std::vector or an std::map. Only when it is an std::vector, it should be resized.
+    // std::map does not have a resize function.
+    points.resize(numberOfPoints);
+  }
+  else
+  {
+    static_assert(std::is_same_v<STLContainerType, std::map<PointIdentifier, PointType>>);
+  }
+
+  auto coordinateIterator = coordinates.cbegin();
+
+  for (PointIdentifier pointIdentifier{}; pointIdentifier < numberOfPoints; ++pointIdentifier)
+  {
+    PointType & point = points[pointIdentifier];
+    std::copy_n(coordinateIterator, PointDimension, point.begin());
+    coordinateIterator += PointDimension;
+  }
+
+  this->Modified();
+}
+
+
 template <typename TPixelType, unsigned int VDimension, typename TMeshTraits>
 auto
 PointSet<TPixelType, VDimension, TMeshTraits>::GetPoints() -> PointsContainer *


Test patch:
diff --git a/Modules/Core/Common/test/CMakeLists.txt b/Modules/Core/Common/test/CMakeLists.txt
index 8f6c1a36334..461b2f96ed7 100644
--- a/Modules/Core/Common/test/CMakeLists.txt
+++ b/Modules/Core/Common/test/CMakeLists.txt
@@ -1748,6 +1748,7 @@ set(ITKCommonGTests
     itkOffsetGTest.cxx
     itkOptimizerParametersGTest.cxx
     itkPointGTest.cxx
+    itkPointSetGTest.cxx
     itkRGBAPixelGTest.cxx
     itkRGBPixelGTest.cxx
     itkShapedImageNeighborhoodRangeGTest.cxx
diff --git a/Modules/Core/Common/test/itkPointSetGTest.cxx b/Modules/Core/Common/test/itkPointSetGTest.cxx
new file mode 100644
index 00000000000..6c1cb8f2ea5
--- /dev/null
+++ b/Modules/Core/Common/test/itkPointSetGTest.cxx
@@ -0,0 +1,82 @@
+/*=========================================================================
+ *
+ *  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.
+ *
+ *=========================================================================*/
+
+// First include the header file to be tested:
+#include "itkPointSet.h"
+#include "../../QuadEdgeMesh/include/itkQuadEdgeMeshTraits.h"
+#include <gtest/gtest.h>
+#include <algorithm> // For equal.
+#include <numeric>   // For iota.
+
+namespace
+{
+template <typename TPointSet>
+void
+TestSetPointsByCoordinates(TPointSet & pointSet)
+{
+  using CoordRepType = typename TPointSet::CoordRepType;
+
+  static constexpr auto PointDimension = TPointSet::PointDimension;
+
+  for (unsigned int numberOfCoordinates{ 1 }; numberOfCoordinates < PointDimension; ++numberOfCoordinates)
+  {
+    // SetPointsByCoordinates is expected to throw an exception when the specified number of coordinates is not a
+    // multiple of PointDimension.
+    EXPECT_THROW(pointSet.SetPointsByCoordinates(std::vector<CoordRepType>(numberOfCoordinates)), itk::ExceptionObject);
+  }
+
+  for (const unsigned int numberOfPoints : { 2, 1, 0 })
+  {
+    std::vector<CoordRepType> coordinates(numberOfPoints * PointDimension);
+
+    // Just make sure that all coordinates have different values, for the purpose of the test.
+    std::iota(coordinates.begin(), coordinates.end(), CoordRepType());
+    {
+      const auto modifiedTime = pointSet.GetMTime();
+      pointSet.SetPointsByCoordinates(coordinates);
+      EXPECT_GT(pointSet.GetMTime(), modifiedTime);
+    }
+
+    using PointsContainerType = typename TPointSet::PointsContainer;
+    using PointIdentifier = typename TPointSet::PointIdentifier;
+    using PointType = typename TPointSet::PointType;
+
+    const typename PointsContainerType::ConstPointer points = pointSet.GetPoints();
+
+    ASSERT_NE(points, nullptr);
+    ASSERT_EQ(points->size(), numberOfPoints);
+
+    const typename PointsContainerType::STLContainerType & stlContainer = points->CastToSTLConstContainer();
+    auto                                                   coordinateIterator = coordinates.cbegin();
+
+    for (PointIdentifier pointIdentifier{}; pointIdentifier < numberOfPoints; ++pointIdentifier)
+    {
+      const PointType & point = stlContainer.at(pointIdentifier);
+      EXPECT_TRUE(std::equal(point.cbegin(), point.cend(), coordinateIterator));
+      coordinateIterator += PointDimension;
+    }
+  }
+}
+} // namespace
+
+
+TEST(PointSet, SetPointsByCoordinates)
+{
+  TestSetPointsByCoordinates(*itk::PointSet<int>::New());
+  TestSetPointsByCoordinates(*itk::PointSet<double, 2, itk::QuadEdgeMeshTraits<double, 2, bool, bool>>::New());
+}