Add methods to BoundaryState to query specific offsets.

The existing functionality worked on radii, but it's helpful to
know if a specific sample location will be in bounds, too.

docs/changelog/BoundaryState-refactor.md

@@ -0,0 +1,145 @@
+Add ability to test exact neighbor offset locations in BoundaryState.
+
+The following methods:
+```
+BoundaryState::InXBoundary
+BoundaryState::InYBoundary
+BoundaryState::InZBoundary
+BoundaryState::InBoundary
+```
+
+have been renamed to:
+
+```
+BoundaryState::IsRadiusInXBoundary
+BoundaryState::IsRadiusInYBoundary
+BoundaryState::IsRadiusInZBoundary
+BoundaryState::IsRadiusInBoundary
+```
+
+to distinguish them from the new methods:
+
+```
+BoundaryState::IsNeighborInXBoundary
+BoundaryState::IsNeighborInYBoundary
+BoundaryState::IsNeighborInZBoundary
+BoundaryState::IsNeighborInBoundary
+```
+
+which check a specific neighbor sample offset instead of a full radius.
+
+The method `BoundaryState::ClampNeighborIndex` has also been added, which clamps
+a 3D neighbor offset vector to the dataset boundaries.
+
+This allows iteration through only the valid points in a neighborhood using
+either of the following patterns:
+
+Using `ClampNeighborIndex` to restrict the iteration space:
+```
+struct MyWorklet : public vtkm::worklet::WorkletPointNeighborhood
+{
+public:
+  using ControlSignature = void(CellSetIn, FieldInNeighborhood, FieldOut);
+  using ExecutionSignature = void(_2, Boundary, _3);
+
+  template <typename InNeighborhoodT, typename OutDataT>
+  VTKM_EXEC void operator()(const InNeighborhoodT& inData,
+                            const vtkm::exec::BoundaryState &boundary,
+                            OutDataT& outData) const
+  {
+    // Clamp the radius to the dataset bounds (discard out-of-bounds points).
+    const auto minRadius = boundary.ClampNeighborIndex({-10, -10, -10});
+    const auto maxRadius = boundary.ClampNeighborIndex({10, 10, 10});
+
+    for (vtkm::IdComponent k = minRadius[2]; k <= maxRadius[2]; ++k)
+    {
+      for (vtkm::IdComponent j = minRadius[1]; j <= maxRadius[1]; ++j)
+      {
+        for (vtkm::IdComponent i = minRadius[0]; i <= maxRadius[0]; ++i)
+        {
+          outData = doSomeConvolution(i, j, k, outdata, inData.Get(i, j, k));
+        }
+      }
+    }
+  }
+};
+```
+
+or, using `IsNeighborInBoundary` methods to skip out-of-bounds loops:
+
+```
+struct MyWorklet : public vtkm::worklet::WorkletPointNeighborhood
+{
+public:
+  using ControlSignature = void(CellSetIn, FieldInNeighborhood, FieldOut);
+  using ExecutionSignature = void(_2, Boundary, _3);
+
+  template <typename InNeighborhoodT, typename OutDataT>
+  VTKM_EXEC void operator()(const InNeighborhoodT& inData,
+                            const vtkm::exec::BoundaryState &boundary,
+                            OutDataT& outData) const
+  {
+    for (vtkm::IdComponent k = -10; k <= 10; ++k)
+    {
+      if (!boundary.IsNeighborInZBoundary(k))
+      {
+        continue;
+      }
+
+      for (vtkm::IdComponent j = -10; j <= 10; ++j)
+      {
+        if (!boundary.IsNeighborInYBoundary(j))
+        {
+          continue;
+        }
+
+        for (vtkm::IdComponent i = -10; i <= 10; ++i)
+        {
+          if (!boundary.IsNeighborInXBoundary(i))
+          {
+            continue;
+          }
+
+          outData = doSomeConvolution(i, j, k, outdata, inData.Get(i, j, k));
+        }
+      }
+    }
+  }
+};
+```
+
+The latter is useful for implementing a convolution that substitutes a constant
+value for out-of-bounds indices:
+
+```
+struct MyWorklet : public vtkm::worklet::WorkletPointNeighborhood
+{
+public:
+  using ControlSignature = void(CellSetIn, FieldInNeighborhood, FieldOut);
+  using ExecutionSignature = void(_2, Boundary, _3);
+
+  template <typename InNeighborhoodT, typename OutDataT>
+  VTKM_EXEC void operator()(const InNeighborhoodT& inData,
+                            const vtkm::exec::BoundaryState &boundary,
+                            OutDataT& outData) const
+  {
+    for (vtkm::IdComponent k = -10; k <= 10; ++k)
+    {
+      for (vtkm::IdComponent j = -10; j <= 10; ++j)
+      {
+        for (vtkm::IdComponent i = -10; i <= 10; ++i)
+        {
+          if (boundary.IsNeighborInBoundary({i, j, k}))
+          {
+            outData = doSomeConvolution(i, j, k, outdata, inData.Get(i, j, k));
+          }
+          else
+          { // substitute zero for out-of-bounds samples:
+            outData = doSomeConvolution(i, j, k, outdata, 0);
+          }
+        }
+      }
+    }
+  }
+};
+```

vtkm/exec/BoundaryState.h

@@ -10,6 +10,7 @@
 #ifndef vtk_m_exec_BoundaryState_h
 #define vtk_m_exec_BoundaryState_h
 
+#include <vtkm/Assert.h>
 #include <vtkm/Math.h>
 
 namespace vtkm
@@ -38,7 +39,7 @@ struct BoundaryState
 
   //@{
   /// Returns true if a neighborhood of the given radius is contained within the bounds of the cell
-  /// set in the X, Y, or Z direction. Returns false if the neighborhood extends ouside of the
+  /// set in the X, Y, or Z direction. Returns false if the neighborhood extends outside of the
   /// boundary of the data in the X, Y, or Z direction.
   ///
   /// The radius defines the size of the neighborhood in terms of how far away it extends from the
@@ -46,22 +47,25 @@ struct BoundaryState
   /// each direction and is 3x3x3. If there is a radius of 2, the neighborhood extends 2 units for
   /// a size of 5x5x5.
   ///
-  VTKM_EXEC bool InXBoundary(vtkm::IdComponent radius) const
+  VTKM_EXEC bool IsRadiusInXBoundary(vtkm::IdComponent radius) const
   {
+    VTKM_ASSERT(radius >= 0);
     return (((this->IJK[0] - radius) >= 0) && ((this->IJK[0] + radius) < this->PointDimensions[0]));
   }
-  VTKM_EXEC bool InYBoundary(vtkm::IdComponent radius) const
+  VTKM_EXEC bool IsRadiusInYBoundary(vtkm::IdComponent radius) const
   {
+    VTKM_ASSERT(radius >= 0);
     return (((this->IJK[1] - radius) >= 0) && ((this->IJK[1] + radius) < this->PointDimensions[1]));
   }
-  VTKM_EXEC bool InZBoundary(vtkm::IdComponent radius) const
+  VTKM_EXEC bool IsRadiusInZBoundary(vtkm::IdComponent radius) const
   {
+    VTKM_ASSERT(radius >= 0);
     return (((this->IJK[2] - radius) >= 0) && ((this->IJK[2] + radius) < this->PointDimensions[2]));
   }
   //@}
 
   /// Returns true if a neighborhood of the given radius is contained within the bounds
-  /// of the cell set. Returns false if the neighborhood extends ouside of the boundary of the
+  /// of the cell set. Returns false if the neighborhood extends outside of the boundary of the
   /// data.
   ///
   /// The radius defines the size of the neighborhood in terms of how far away it extends from the
@@ -69,15 +73,47 @@ struct BoundaryState
   /// each direction and is 3x3x3. If there is a radius of 2, the neighborhood extends 2 units for
   /// a size of 5x5x5.
   ///
-  VTKM_EXEC bool InBoundary(vtkm::IdComponent radius) const
+  VTKM_EXEC bool IsRadiusInBoundary(vtkm::IdComponent radius) const
   {
-    return this->InXBoundary(radius) && this->InYBoundary(radius) && this->InZBoundary(radius);
+    return this->IsRadiusInXBoundary(radius) && this->IsRadiusInYBoundary(radius) &&
+      this->IsRadiusInZBoundary(radius);
+  }
+
+  //@{
+  /// Returns true if the neighbor at the specified @a offset is contained
+  /// within the bounds of the cell set in the X, Y, or Z direction. Returns
+  /// false if the neighbor falls outside of the boundary of the data in the X,
+  /// Y, or Z direction.
+  ///
+  VTKM_EXEC bool IsNeighborInXBoundary(vtkm::IdComponent offset) const
+  {
+    return (((this->IJK[0] + offset) >= 0) && ((this->IJK[0] + offset) < this->PointDimensions[0]));
+  }
+  VTKM_EXEC bool IsNeighborInYBoundary(vtkm::IdComponent offset) const
+  {
+    return (((this->IJK[1] + offset) >= 0) && ((this->IJK[1] + offset) < this->PointDimensions[1]));
+  }
+  VTKM_EXEC bool IsNeighborInZBoundary(vtkm::IdComponent offset) const
+  {
+    return (((this->IJK[2] + offset) >= 0) && ((this->IJK[2] + offset) < this->PointDimensions[2]));
+  }
+  //@}
+
+  /// Returns true if the neighbor at the specified offset vector is contained
+  /// within the bounds of the cell set. Returns false if the neighbor falls
+  /// outside of the boundary of the data.
+  ///
+  VTKM_EXEC bool IsNeighborInBoundary(const vtkm::Vec<vtkm::IdComponent, 3>& neighbor) const
+  {
+    return this->IsNeighborInXBoundary(neighbor[0]) && this->IsNeighborInYBoundary(neighbor[1]) &&
+      this->IsNeighborInZBoundary(neighbor[2]);
   }
 
   /// Returns the minimum neighborhood indices that are within the bounds of the data.
   ///
   VTKM_EXEC vtkm::Vec<vtkm::IdComponent, 3> MinNeighborIndices(vtkm::IdComponent radius) const
   {
+    VTKM_ASSERT(radius >= 0);
     vtkm::Vec<vtkm::IdComponent, 3> minIndices;
 
     for (vtkm::IdComponent component = 0; component < 3; ++component)
@@ -99,6 +135,7 @@ struct BoundaryState
   ///
   VTKM_EXEC vtkm::Vec<vtkm::IdComponent, 3> MaxNeighborIndices(vtkm::IdComponent radius) const
   {
+    VTKM_ASSERT(radius >= 0);
     vtkm::Vec<vtkm::IdComponent, 3> maxIndices;
 
     for (vtkm::IdComponent component = 0; component < 3; ++component)
@@ -117,9 +154,6 @@ struct BoundaryState
     return maxIndices;
   }
 
-  //todo: This needs to work with BoundaryConstantValue
-  //todo: This needs to work with BoundaryPeroidic
-
   //@{
   /// Takes a local neighborhood index (in the ranges of -neighborhood size to neighborhood size)
   /// and returns the ijk of the equivalent point in the full data set. If the given value is out
@@ -143,8 +177,29 @@ struct BoundaryState
   }
   //@}
 
-  //todo: This needs to work with BoundaryConstantValue
+  //@{
-  //todo: This needs to work with BoundaryPeroidic
+  /// Takes a local neighborhood index (in the ranges of -neighborhood size to
+  /// neighborhood size), clamps it to the dataset bounds, and returns a new
+  /// neighborhood index. For example, if given a neighbor index that is past
+  /// the minimum x range of the data, the neighbor index of the minimum x
+  /// boundary is returned.
+  ///
+  VTKM_EXEC vtkm::Vec<vtkm::IdComponent, 3> ClampNeighborIndex(
+    const vtkm::Vec<vtkm::IdComponent, 3>& neighbor) const
+  {
+    const vtkm::Id3 fullIndex = this->IJK + neighbor;
+    const vtkm::Id3 clampedFullIndex =
+      vtkm::Max(vtkm::Id3(0), vtkm::Min(this->PointDimensions - vtkm::Id3(1), fullIndex));
+    return vtkm::Vec<vtkm::IdComponent, 3>{ clampedFullIndex - this->IJK };
+  }
+
+  VTKM_EXEC vtkm::Vec<vtkm::IdComponent, 3> ClampNeighborIndex(vtkm::IdComponent neighborI,
+                                                               vtkm::IdComponent neighborJ,
+                                                               vtkm::IdComponent neighborK) const
+  {
+    return this->ClampNeighborIndex(vtkm::make_Vec(neighborI, neighborJ, neighborK));
+  }
+  //@}
 
   //@{
   /// Takes a local neighborhood index (in the ranges of -neighborhood size to neighborhood size)

vtkm/exec/arg/testing/UnitTestFetchArrayNeighborhoodIn.cxx

@@ -44,23 +44,37 @@ void verify_neighbors(NeighborhoodType neighbors, vtkm::Id index, vtkm::Id3 inde
 
   //Verify the boundary flags first
   VTKM_TEST_ASSERT(((index3d[0] != 0) && (index3d[0] != (POINT_DIMS[0] - 1))) ==
-                     boundary->InXBoundary(1),
+                     boundary->IsRadiusInXBoundary(1),
-                   "Got invalid X boundary");
+                   "Got invalid X radius boundary");
   VTKM_TEST_ASSERT(((index3d[1] != 0) && (index3d[1] != (POINT_DIMS[1] - 1))) ==
-                     boundary->InYBoundary(1),
+                     boundary->IsRadiusInYBoundary(1),
-                   "Got invalid Y boundary");
+                   "Got invalid Y radius boundary");
   VTKM_TEST_ASSERT(((index3d[2] != 0) && (index3d[2] != (POINT_DIMS[2] - 1))) ==
-                     boundary->InZBoundary(1),
+                     boundary->IsRadiusInZBoundary(1),
-                   "Got invalid Z boundary");
+                   "Got invalid Z radius boundary");
+
+  VTKM_TEST_ASSERT((index3d[0] != 0) == boundary->IsNeighborInXBoundary(-1),
+                   "Got invalid X negative neighbor boundary");
+  VTKM_TEST_ASSERT((index3d[1] != 0) == boundary->IsNeighborInYBoundary(-1),
+                   "Got invalid Y negative neighbor boundary");
+  VTKM_TEST_ASSERT((index3d[2] != 0) == boundary->IsNeighborInZBoundary(-1),
+                   "Got invalid Z negative neighbor boundary");
+
+  VTKM_TEST_ASSERT((index3d[0] != (POINT_DIMS[0] - 1)) == boundary->IsNeighborInXBoundary(1),
+                   "Got invalid X positive neighbor boundary");
+  VTKM_TEST_ASSERT((index3d[1] != (POINT_DIMS[1] - 1)) == boundary->IsNeighborInYBoundary(1),
+                   "Got invalid Y positive neighbor boundary");
+  VTKM_TEST_ASSERT((index3d[2] != (POINT_DIMS[2] - 1)) == boundary->IsNeighborInZBoundary(1),
+                   "Got invalid Z positive neighbor boundary");
 
   VTKM_TEST_ASSERT(((boundary->MinNeighborIndices(1)[0] == -1) &&
-                    (boundary->MaxNeighborIndices(1)[0] == 1)) == boundary->InXBoundary(1),
+                    (boundary->MaxNeighborIndices(1)[0] == 1)) == boundary->IsRadiusInXBoundary(1),
                    "Got invalid min/max X indices");
   VTKM_TEST_ASSERT(((boundary->MinNeighborIndices(1)[1] == -1) &&
-                    (boundary->MaxNeighborIndices(1)[1] == 1)) == boundary->InYBoundary(1),
+                    (boundary->MaxNeighborIndices(1)[1] == 1)) == boundary->IsRadiusInYBoundary(1),
                    "Got invalid min/max Y indices");
   VTKM_TEST_ASSERT(((boundary->MinNeighborIndices(1)[2] == -1) &&
-                    (boundary->MaxNeighborIndices(1)[2] == 1)) == boundary->InZBoundary(1),
+                    (boundary->MaxNeighborIndices(1)[2] == 1)) == boundary->IsRadiusInZBoundary(1),
                    "Got invalid min/max Z indices");
 
   T forwardX = neighbors.Get(1, 0, 0);

vtkm/worklet/gradient/StructuredPointGradient.h

@@ -57,9 +57,9 @@ struct StructuredPointGradient : public vtkm::worklet::WorkletPointNeighborhood
     T deta = inputField.Get(0, 1, 0) - inputField.Get(0, -1, 0);
     T dzeta = inputField.Get(0, 0, 1) - inputField.Get(0, 0, -1);
 
-    dxi = (boundary.InXBoundary(1) ? dxi * 0.5f : dxi);
+    dxi = (boundary.IsRadiusInXBoundary(1) ? dxi * 0.5f : dxi);
-    deta = (boundary.InYBoundary(1) ? deta * 0.5f : deta);
+    deta = (boundary.IsRadiusInYBoundary(1) ? deta * 0.5f : deta);
-    dzeta = (boundary.InZBoundary(1) ? dzeta * 0.5f : dzeta);
+    dzeta = (boundary.IsRadiusInZBoundary(1) ? dzeta * 0.5f : dzeta);
 
     outputGradient[0] = static_cast<OT>(xi[0] * dxi + eta[0] * deta + zeta[0] * dzeta);
     outputGradient[1] = static_cast<OT>(xi[1] * dxi + eta[1] * deta + zeta[1] * dzeta);
@@ -84,9 +84,9 @@ struct StructuredPointGradient : public vtkm::worklet::WorkletPointNeighborhood
 
     CoordType r = inputPoints.Portal.GetSpacing();
 
-    r[0] = (boundary.InXBoundary(1) ? r[0] * 0.5f : r[0]);
+    r[0] = (boundary.IsRadiusInXBoundary(1) ? r[0] * 0.5f : r[0]);
-    r[1] = (boundary.InYBoundary(1) ? r[1] * 0.5f : r[1]);
+    r[1] = (boundary.IsRadiusInYBoundary(1) ? r[1] * 0.5f : r[1]);
-    r[2] = (boundary.InZBoundary(1) ? r[2] * 0.5f : r[2]);
+    r[2] = (boundary.IsRadiusInZBoundary(1) ? r[2] * 0.5f : r[2]);
 
     const T dx = inputField.Get(1, 0, 0) - inputField.Get(-1, 0, 0);
     const T dy = inputField.Get(0, 1, 0) - inputField.Get(0, -1, 0);
@@ -113,9 +113,9 @@ struct StructuredPointGradient : public vtkm::worklet::WorkletPointNeighborhood
     CoordType eta = inputPoints.Get(0, 1, 0) - inputPoints.Get(0, -1, 0);
     CoordType zeta = inputPoints.Get(0, 0, 1) - inputPoints.Get(0, 0, -1);
 
-    xi = (boundary.InXBoundary(1) ? xi * 0.5f : xi);
+    xi = (boundary.IsRadiusInXBoundary(1) ? xi * 0.5f : xi);
-    eta = (boundary.InYBoundary(1) ? eta * 0.5f : eta);
+    eta = (boundary.IsRadiusInYBoundary(1) ? eta * 0.5f : eta);
-    zeta = (boundary.InZBoundary(1) ? zeta * 0.5f : zeta);
+    zeta = (boundary.IsRadiusInZBoundary(1) ? zeta * 0.5f : zeta);
 
     CT aj = xi[0] * eta[1] * zeta[2] + xi[1] * eta[2] * zeta[0] + xi[2] * eta[0] * zeta[1] -
       xi[2] * eta[1] * zeta[0] - xi[1] * eta[0] * zeta[2] - xi[0] * eta[2] * zeta[1];

vtkm/worklet/testing/UnitTestWorkletMapPointNeighborhood.cxx

@@ -44,26 +44,50 @@ struct MaxNeighborValue : public vtkm::worklet::WorkletPointNeighborhood
 
     auto* nboundary = inputField.Boundary;
 
-    if (!(nboundary->InXBoundary(1) == boundary.InXBoundary(1)))
+    if (!(nboundary->IsRadiusInXBoundary(1) == boundary.IsRadiusInXBoundary(1)))
     {
       this->RaiseError("Got invalid XPos boundary state");
     }
 
-    if (!(nboundary->InYBoundary(1) == boundary.InYBoundary(1)))
+    if (!(nboundary->IsRadiusInYBoundary(1) == boundary.IsRadiusInYBoundary(1)))
     {
       this->RaiseError("Got invalid YPos boundary state");
     }
 
-    if (!(nboundary->InZBoundary(1) == boundary.InZBoundary(1)))
+    if (!(nboundary->IsRadiusInZBoundary(1) == boundary.IsRadiusInZBoundary(1)))
     {
       this->RaiseError("Got invalid ZPos boundary state");
     }
 
-    if (!(nboundary->InBoundary(1) == boundary.InBoundary(1)))
+    if (!(nboundary->IsRadiusInBoundary(1) == boundary.IsRadiusInBoundary(1)))
     {
       this->RaiseError("Got invalid boundary state");
     }
 
+    if (nboundary->IsRadiusInXBoundary(1) !=
+        (boundary.IsNeighborInXBoundary(-1) && boundary.IsNeighborInXBoundary(1)))
+    {
+      this->RaiseError("Neighbor/Radius boundary mismatch in X dimension.");
+    }
+
+    if (nboundary->IsRadiusInYBoundary(1) !=
+        (boundary.IsNeighborInYBoundary(-1) && boundary.IsNeighborInYBoundary(1)))
+    {
+      this->RaiseError("Neighbor/Radius boundary mismatch in Y dimension.");
+    }
+
+    if (nboundary->IsRadiusInZBoundary(1) !=
+        (boundary.IsNeighborInZBoundary(-1) && boundary.IsNeighborInZBoundary(1)))
+    {
+      this->RaiseError("Neighbor/Radius boundary mismatch in Z dimension.");
+    }
+
+    if (nboundary->IsRadiusInBoundary(1) !=
+        (boundary.IsNeighborInBoundary({ -1 }) && boundary.IsNeighborInBoundary({ 1 })))
+    {
+      this->RaiseError("Neighbor/Radius boundary mismatch.");
+    }
+
 
     auto minNeighbors = boundary.MinNeighborIndices(1);
     auto maxNeighbors = boundary.MaxNeighborIndices(1);