Merge topic 'density-interp-any-field'

0ba818e90 Enable any scalar field in ParticleDensity filters
d38efa87f Allow fields of any type in NDHistogram

Acked-by: Kitware Robot <kwrobot@kitware.com>
Acked-by: Li-Ta Lo <ollie@lanl.gov>
Merge-request: !2978

docs/changelog/interp-any-field.md

@@ -1,4 +1,4 @@
-# Update filters' field map to work on any field type
+# Update filters' field map and execute to work on any field type
 
 Several filters implemented their map field by checking for common field
 types and interpolated those. Although there was a float fallback to catch
@@ -10,7 +10,8 @@ so it is still possible to drop fields.
 The map field functions for these filters have been changed to support all
 possible types. This is done by using the extract component functionality
 to get data from any type of array. The following filters have been
-updated.
+updated. In some circumstances where it makes sense, a simple float
+fallback is used.
 
   * `CleanGrid`
   * `CellAverage`
@@ -19,3 +20,6 @@ updated.
   * `Contour`
   * `MIRFilter`
   * `PointAverage`
+  * `NDHistogram`
+  * `ParticleDensityCloudInCell`
+  * `ParticleDensityNearestGridPoint`

vtkm/filter/density_estimate/ParticleDensityCloudInCell.cxx

@@ -111,28 +111,15 @@ VTKM_CONT vtkm::cont::DataSet ParticleDensityCloudInCell::DoExecute(const cont::
     uniform.AddField(vtkm::cont::make_FieldPoint("density", density));
   };
 
-  // Note: This is the so called Immediately-Invoked Function Expression (IIFE). Here we define
-  // a lambda expression and immediately call it at the end. This allows us to not declare an
-  // UnknownArrayHandle first and then assign it in the if-else statement. If I really want to
-  // show-off, I can even inline the `fieldArray` variable and turn it into a long expression.
-  auto fieldArray = [&]() -> vtkm::cont::UnknownArrayHandle {
-    if (this->ComputeNumberDensity)
-    {
-      return vtkm::cont::make_ArrayHandleConstant(vtkm::FloatDefault{ 1 },
-                                                  input.GetNumberOfPoints());
-    }
-    else
-    {
-      return this->GetFieldFromDataSet(input).GetData();
-    }
-  }();
-  fieldArray.CastAndCallForTypes<
-    vtkm::TypeListFieldScalar,
-    vtkm::ListAppend<VTKM_DEFAULT_STORAGE_LIST, vtkm::List<vtkm::cont::StorageTagConstant>>>(
-    resolveType);
-
-  // Deposition of the input field to the output field is already mapping. No need to map other
-  // fields.
+  if (this->ComputeNumberDensity)
+  {
+    resolveType(
+      vtkm::cont::make_ArrayHandleConstant(vtkm::FloatDefault{ 1 }, input.GetNumberOfPoints()));
+  }
+  else
+  {
+    this->CastAndCallScalarField(this->GetFieldFromDataSet(input), resolveType);
+  }
   return uniform;
 }
 } // namespace density_estimate

vtkm/filter/density_estimate/ParticleDensityNearestGridPoint.cxx

@@ -95,25 +95,15 @@ VTKM_CONT vtkm::cont::DataSet ParticleDensityNearestGridPoint::DoExecute(
     uniform.AddField(vtkm::cont::make_FieldCell("density", density));
   };
 
-  // Note: This is the so called Immediately-Invoked Function Expression (IIFE). Here we define
-  // a lambda expression and immediately call it at the end. This allows us to not declare an
-  // UnknownArrayHandle first and then assign it in the if-else statement. If I really want to
-  // show-off, I can even inline the `fieldArray` variable and turn it into a long expression.
-  auto fieldArray = [&]() -> vtkm::cont::UnknownArrayHandle {
-    if (this->ComputeNumberDensity)
-    {
-      return vtkm::cont::make_ArrayHandleConstant(vtkm::FloatDefault{ 1 },
-                                                  input.GetNumberOfPoints());
-    }
-    else
-    {
-      return this->GetFieldFromDataSet(input).GetData();
-    }
-  }();
-  fieldArray.CastAndCallForTypes<
-    vtkm::TypeListFieldScalar,
-    vtkm::ListAppend<VTKM_DEFAULT_STORAGE_LIST, vtkm::List<vtkm::cont::StorageTagConstant>>>(
-    resolveType);
+  if (this->ComputeNumberDensity)
+  {
+    resolveType(
+      vtkm::cont::make_ArrayHandleConstant(vtkm::FloatDefault{ 1 }, input.GetNumberOfPoints()));
+  }
+  else
+  {
+    this->CastAndCallScalarField(this->GetFieldFromDataSet(input), resolveType);
+  }
 
   // Deposition of the input field to the output field is already mapping. No need to map other
   // fields.

vtkm/filter/density_estimate/worklet/NDimsHistogram.h

@@ -41,29 +41,6 @@ public:
     vtkm::cont::ArrayCopy(constant0Array, Bin1DIndex);
   }
 
-  // Add a field and the bin number for this field
-  // Return: rangeOfRange is min max value of this array
-  //         binDelta is delta of a bin
-  template <typename HandleType>
-  void AddField(const HandleType& fieldArray,
-                vtkm::Id numberOfBins,
-                vtkm::Range& rangeOfValues,
-                vtkm::Float64& binDelta)
-  {
-    NumberOfBins.push_back(numberOfBins);
-
-    if (fieldArray.GetNumberOfValues() != NumDataPoints)
-    {
-      throw vtkm::cont::ErrorBadValue("Array lengths does not match");
-    }
-    else
-    {
-      vtkm::cont::CastAndCall(
-        fieldArray.ResetTypes(vtkm::TypeListScalarAll{}, VTKM_DEFAULT_STORAGE_LIST{}),
-        vtkm::worklet::histogram::ComputeBins(Bin1DIndex, numberOfBins, rangeOfValues, binDelta));
-    }
-  }
-
   // Add a field and the bin number for this field along with specific range of the data
   // Return: binDelta is delta of a bin
   template <typename HandleType>
@@ -71,7 +48,7 @@ public:
                 vtkm::Id numberOfBins,
                 vtkm::Range& rangeOfValues,
                 vtkm::Float64& binDelta,
-                bool rangeProvided)
+                bool rangeProvided = false)
   {
     NumberOfBins.push_back(numberOfBins);
 
@@ -81,9 +58,35 @@ public:
     }
     else
     {
-      CastAndCall(fieldArray.ResetTypes(vtkm::TypeListScalarAll()),
-                  vtkm::worklet::histogram::ComputeBins(
-                    Bin1DIndex, numberOfBins, rangeOfValues, binDelta, rangeProvided));
+      auto computeBins = [&](auto resolvedField) {
+        // CastAndCallWithExtractedArray should give us an ArrayHandleRecombineVec
+        using T = typename std::decay_t<decltype(resolvedField)>::ValueType::ComponentType;
+        vtkm::cont::ArrayHandleRecombineVec<T> recombineField{ resolvedField };
+        if (recombineField.GetNumberOfComponents() != 1)
+        {
+          VTKM_LOG_S(vtkm::cont::LogLevel::Warn,
+                     "NDHistogram expects scalar fields, but was given field with "
+                       << recombineField.GetNumberOfComponents()
+                       << " components. Extracting first component.");
+        }
+        vtkm::cont::ArrayHandleStride<T> field =
+          vtkm::cont::ArrayExtractComponent(recombineField, 0);
+        if (!rangeProvided)
+        {
+          const vtkm::Vec<T, 2> initValue(vtkm::cont::ArrayGetValue(0, field));
+          vtkm::Vec<T, 2> minMax =
+            vtkm::cont::Algorithm::Reduce(field, initValue, vtkm::MinAndMax<T>());
+          rangeOfValues.Min = static_cast<vtkm::Float64>(minMax[0]);
+          rangeOfValues.Max = static_cast<vtkm::Float64>(minMax[1]);
+        }
+        binDelta = vtkm::worklet::histogram::compute_delta(
+          rangeOfValues.Min, rangeOfValues.Max, numberOfBins);
+
+        vtkm::worklet::histogram::SetHistogramBin<T> binWorklet(
+          numberOfBins, rangeOfValues.Min, binDelta);
+        vtkm::cont::Invoker{}(binWorklet, field, this->Bin1DIndex, this->Bin1DIndex);
+      };
+      fieldArray.CastAndCallWithExtractedArray(computeBins);
     }
   }
 

vtkm/filter/density_estimate/worklet/histogram/ComputeNDHistogram.h

@@ -78,63 +78,6 @@ public:
   }
 };
 
-class ComputeBins
-{
-public:
-  VTKM_CONT
-  ComputeBins(vtkm::cont::ArrayHandle<vtkm::Id>& _bin1DIdx,
-              vtkm::Id& _numOfBins,
-              vtkm::Range& _minMax,
-              vtkm::Float64& _binDelta)
-    : Bin1DIdx(_bin1DIdx)
-    , NumOfBins(_numOfBins)
-    , MinMax(_minMax)
-    , BinDelta(_binDelta)
-    , RangeProvided(false)
-  {
-  }
-
-  VTKM_CONT
-  ComputeBins(vtkm::cont::ArrayHandle<vtkm::Id>& _bin1DIdx,
-              vtkm::Id& _numOfBins,
-              vtkm::Range& _minMax,
-              vtkm::Float64& _binDelta,
-              bool _rangeProvided)
-    : Bin1DIdx(_bin1DIdx)
-    , NumOfBins(_numOfBins)
-    , MinMax(_minMax)
-    , BinDelta(_binDelta)
-    , RangeProvided(_rangeProvided)
-  {
-  }
-
-  template <typename T, typename Storage>
-  VTKM_CONT void operator()(const vtkm::cont::ArrayHandle<T, Storage>& field) const
-  {
-    if (!RangeProvided)
-    {
-      const vtkm::Vec<T, 2> initValue(vtkm::cont::ArrayGetValue(0, field));
-      vtkm::Vec<T, 2> minMax =
-        vtkm::cont::Algorithm::Reduce(field, initValue, vtkm::MinAndMax<T>());
-      MinMax.Min = static_cast<vtkm::Float64>(minMax[0]);
-      MinMax.Max = static_cast<vtkm::Float64>(minMax[1]);
-    }
-    BinDelta = compute_delta(MinMax.Min, MinMax.Max, NumOfBins);
-
-    SetHistogramBin<T> binWorklet(NumOfBins, MinMax.Min, BinDelta);
-    vtkm::worklet::DispatcherMapField<vtkm::worklet::histogram::SetHistogramBin<T>>
-      setHistogramBinDispatcher(binWorklet);
-    setHistogramBinDispatcher.Invoke(field, Bin1DIdx, Bin1DIdx);
-  }
-
-private:
-  vtkm::cont::ArrayHandle<vtkm::Id>& Bin1DIdx;
-  vtkm::Id& NumOfBins;
-  vtkm::Range& MinMax;
-  vtkm::Float64& BinDelta;
-  bool RangeProvided;
-};
-
 // Convert N-dims bin index into 1D index
 class ConvertHistBinToND : public vtkm::worklet::WorkletMapField
 {