diff --git a/benchmarking/BenchmarkFilters.cxx b/benchmarking/BenchmarkFilters.cxx index d0d37f8ee..dfe03e3bb 100644 --- a/benchmarking/BenchmarkFilters.cxx +++ b/benchmarking/BenchmarkFilters.cxx @@ -46,8 +46,8 @@ #include +#include #include -#include #include #include // for std::tolower @@ -1357,10 +1357,10 @@ int BenchmarkBody(int argc, char** argv, const vtkm::cont::InitializeResult& con { std::cout << "Generating " << waveletDim << "x" << waveletDim << "x" << waveletDim << " wavelet...\n"; - vtkm::worklet::WaveletGenerator gen; - gen.SetExtent({ 0 }, { waveletDim }); + vtkm::source::Wavelet source; + source.SetExtent({ 0 }, { waveletDim }); - InputDataSet = gen.GenerateDataSet(config.Device); + InputDataSet = source.Execute(); } if (tetra) diff --git a/benchmarking/CMakeLists.txt b/benchmarking/CMakeLists.txt index cc95b4da3..201492b01 100644 --- a/benchmarking/CMakeLists.txt +++ b/benchmarking/CMakeLists.txt @@ -42,7 +42,7 @@ set(benchmarks ) foreach (benchmark ${benchmarks}) - add_benchmark(NAME ${benchmark} FILE ${benchmark}.cxx LIBS vtkm_filter) + add_benchmark(NAME ${benchmark} FILE ${benchmark}.cxx LIBS vtkm_source vtkm_filter) endforeach () if(TARGET vtkm_rendering) diff --git a/vtkm/filter/testing/UnitTestSplitSharpEdgesFilter.cxx b/vtkm/filter/testing/UnitTestSplitSharpEdgesFilter.cxx index 1d842c2c6..601ea1219 100644 --- a/vtkm/filter/testing/UnitTestSplitSharpEdgesFilter.cxx +++ b/vtkm/filter/testing/UnitTestSplitSharpEdgesFilter.cxx @@ -15,7 +15,7 @@ #include #include -#include +#include namespace { @@ -115,13 +115,11 @@ vtkm::cont::DataSet Make3DExplicitSimpleCube() vtkm::cont::DataSet Make3DWavelet() { - vtkm::worklet::WaveletGenerator wavelet; - wavelet.SetMinimumExtent({ -25 }); - wavelet.SetMaximumExtent({ 25 }); + vtkm::source::Wavelet wavelet({ -25 }, { 25 }); wavelet.SetFrequency({ 60, 30, 40 }); wavelet.SetMagnitude({ 5 }); - vtkm::cont::DataSet result = wavelet.GenerateDataSet(); + vtkm::cont::DataSet result = wavelet.Execute(); return result; } diff --git a/vtkm/source/CMakeLists.txt b/vtkm/source/CMakeLists.txt index 07e78e47b..c9ab20243 100644 --- a/vtkm/source/CMakeLists.txt +++ b/vtkm/source/CMakeLists.txt @@ -11,11 +11,13 @@ set(headers Source.h Tangle.h + Wavelet.h ) set(device_sources Source.cxx Tangle.cxx + Wavelet.cxx ) vtkm_library(NAME vtkm_source diff --git a/vtkm/source/Wavelet.cxx b/vtkm/source/Wavelet.cxx new file mode 100644 index 000000000..8d00ccf22 --- /dev/null +++ b/vtkm/source/Wavelet.cxx @@ -0,0 +1,183 @@ +//============================================================================ +// Copyright (c) Kitware, Inc. +// All rights reserved. +// See LICENSE.txt for details. +// +// This software is distributed WITHOUT ANY WARRANTY; without even +// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR +// PURPOSE. See the above copyright notice for more information. +//============================================================================ + +#include +#include + +namespace +{ +inline vtkm::FloatDefault computeScaleFactor(vtkm::Id min, vtkm::Id max) +{ + return (min < max) ? (1.f / static_cast(max - min)) + : static_cast(1.); +} +} +namespace vtkm +{ +namespace source +{ +namespace wavelet +{ +template +struct Worker : public vtkm::exec::FunctorBase +{ + using OutputHandleType = vtkm::cont::ArrayHandle; + using OutputPortalType = decltype(std::declval().PrepareForOutput(0, Device())); + using Vec3F = vtkm::Vec3f; + + Vec3F Center; + Vec3F Spacing; + Vec3F Frequency; + Vec3F Magnitude; + Vec3F MinimumPoint; + Vec3F Scale; + vtkm::Id3 Offset; + vtkm::Id3 Dims; + vtkm::FloatDefault MaximumValue; + vtkm::FloatDefault Temp2; + OutputPortalType Portal; + + VTKM_CONT + Worker(OutputHandleType& output, + const Vec3F& center, + const Vec3F& spacing, + const Vec3F& frequency, + const Vec3F& magnitude, + const Vec3F& minimumPoint, + const Vec3F& scale, + const vtkm::Id3& offset, + const vtkm::Id3& dims, + vtkm::FloatDefault maximumValue, + vtkm::FloatDefault temp2) + : Center(center) + , Spacing(spacing) + , Frequency(frequency) + , Magnitude(magnitude) + , MinimumPoint(minimumPoint) + , Scale(scale) + , Offset(offset) + , Dims(dims) + , MaximumValue(maximumValue) + , Temp2(temp2) + , Portal(output.PrepareForOutput((dims[0] * dims[1] * dims[2]), Device{})) + { + } + + VTKM_EXEC + void operator()(const vtkm::Id3& ijk) const + { + // map ijk to the point location, accounting for spacing: + const Vec3F loc = Vec3F(ijk + this->Offset) * this->Spacing; + + // Compute the distance from the center of the gaussian: + const Vec3F scaledLoc = (this->Center - loc) * this->Scale; + vtkm::FloatDefault gaussSum = vtkm::Dot(scaledLoc, scaledLoc); + + const Vec3F periodicContribs{ + this->Magnitude[0] * vtkm::Sin(this->Frequency[0] * scaledLoc[0]), + this->Magnitude[1] * vtkm::Sin(this->Frequency[1] * scaledLoc[1]), + this->Magnitude[2] * vtkm::Cos(this->Frequency[2] * scaledLoc[2]), + }; + + // The vtkRTAnalyticSource documentation says the periodic contributions + // should be multiplied in, but the implementation adds them. We'll do as + // they do, not as they say. + const vtkm::FloatDefault scalar = this->MaximumValue * vtkm::Exp(-gaussSum * this->Temp2) + + periodicContribs[0] + periodicContribs[1] + periodicContribs[2]; + + // Compute output location + // (see ConnectivityStructuredInternals<3>::LogicalToFlatPointIndex) + const vtkm::Id scalarIdx = ijk[0] + this->Dims[0] * (ijk[1] + this->Dims[1] * ijk[2]); + this->Portal.Set(scalarIdx, scalar); + } +}; + +struct runWorker +{ + template + inline bool operator()(Device, const vtkm::Id3 dims, Args... args) const + { + using Algo = vtkm::cont::DeviceAdapterAlgorithm; + Worker worker{ args... }; + Algo::Schedule(worker, dims); + return true; + } +}; +} // namespace wavelet + +Wavelet::Wavelet(vtkm::Id3 minExtent, vtkm::Id3 maxExtent) + : Center{ minExtent - ((minExtent - maxExtent) / 2) } + , Spacing{ 1. } + , Frequency{ 60., 30., 40. } + , Magnitude{ 10., 18., 5. } + , MinimumExtent{ minExtent } + , MaximumExtent{ maxExtent } + , MaximumValue{ 255. } + , StandardDeviation{ 0.5 } +{ +} + +vtkm::cont::DataSet Wavelet::Execute() const +{ + VTKM_LOG_SCOPE_FUNCTION(vtkm::cont::LogLevel::Perf); + + // Create points: + const vtkm::Id3 dims{ this->MaximumExtent - this->MinimumExtent + vtkm::Id3{ 1 } }; + const vtkm::Vec3f origin{ this->MinimumExtent }; + vtkm::cont::CoordinateSystem coords{ "coordinates", dims, origin, this->Spacing }; + + // And cells: + vtkm::cont::CellSetStructured<3> cellSet; + cellSet.SetPointDimensions(dims); + + // Compile the dataset: + vtkm::cont::DataSet dataSet; + dataSet.AddCoordinateSystem(coords); + dataSet.SetCellSet(cellSet); + + // Scalars, too + vtkm::cont::Field field = this->GeneratePointField("scalars"); + dataSet.AddField(field); + + return dataSet; +} + +vtkm::cont::Field Wavelet::GeneratePointField(const std::string& name) const +{ + VTKM_LOG_SCOPE_FUNCTION(vtkm::cont::LogLevel::Perf); + + const vtkm::Id3 dims{ this->MaximumExtent - this->MinimumExtent + vtkm::Id3{ 1 } }; + vtkm::Vec3f minPt = vtkm::Vec3f(this->MinimumExtent) * this->Spacing; + vtkm::FloatDefault temp2 = 1.f / (2.f * this->StandardDeviation * this->StandardDeviation); + vtkm::Vec3f scale{ computeScaleFactor(this->MinimumExtent[0], this->MaximumExtent[0]), + computeScaleFactor(this->MinimumExtent[1], this->MaximumExtent[1]), + computeScaleFactor(this->MinimumExtent[2], this->MaximumExtent[2]) }; + + + vtkm::cont::ArrayHandle output; + vtkm::cont::TryExecuteOnDevice(this->Invoke.GetDevice(), + wavelet::runWorker{}, + dims, + output, + this->Center, + this->Spacing, + this->Frequency, + this->Magnitude, + minPt, + scale, + this->MinimumExtent, + dims, + this->MaximumValue, + temp2); + return vtkm::cont::make_FieldPoint(name, output); +} + +} // namespace source +} // namespace vtkm diff --git a/vtkm/source/Wavelet.h b/vtkm/source/Wavelet.h new file mode 100644 index 000000000..bfc6fff41 --- /dev/null +++ b/vtkm/source/Wavelet.h @@ -0,0 +1,111 @@ +//============================================================================ +// Copyright (c) Kitware, Inc. +// All rights reserved. +// See LICENSE.txt for details. +// +// This software is distributed WITHOUT ANY WARRANTY; without even +// the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR +// PURPOSE. See the above copyright notice for more information. +//============================================================================ + +#ifndef vtk_m_source_Wavelet_h +#define vtk_m_source_Wavelet_h + +#include + +namespace vtkm +{ +namespace source +{ +/** + * @brief The Wavelet source creates a dataset similar to VTK's + * vtkRTAnalyticSource. + * + * This class generates a predictable structured dataset with a smooth yet + * interesting set of scalars, which is useful for testing and benchmarking. + * + * The Execute method creates a complete structured dataset that have a + * point field names 'scalars' + * + * The scalars are computed as: + * + * ``` + * MaxVal * Gauss + MagX * sin(FrqX*x) + MagY * sin(FrqY*y) + MagZ * cos(FrqZ*z) + * ``` + * + * The dataset properties are determined by: + * - `Minimum/MaximumExtent`: The logical point extents of the dataset. + * - `Spacing`: The distance between points of the dataset. + * - `Center`: The center of the dataset. + * + * The scalar functions is control via: + * - `Center`: The center of a Gaussian contribution to the scalars. + * - `StandardDeviation`: The unscaled width of a Gaussian contribution. + * - `MaximumValue`: Upper limit of the scalar range. + * - `Frequency`: The Frq[XYZ] parameters of the periodic contributions. + * - `Magnitude`: The Mag[XYZ] parameters of the periodic contributions. + * + * By default, the following parameters are used: + * - `Extents`: { -10, -10, -10 } `-->` { 10, 10, 10 } + * - `Spacing`: { 1, 1, 1 } + * - `Center`: { 0, 0, 0 } + * - `StandardDeviation`: 0.5 + * - `MaximumValue`: 255 + * - `Frequency`: { 60, 30, 40 } + * - `Magnitude`: { 10, 18, 5 } + */ +class VTKM_SOURCE_EXPORT Wavelet final : public vtkm::source::Source +{ +public: + VTKM_CONT + Wavelet(vtkm::Id3 minExtent = { -10 }, vtkm::Id3 maxExtent = { 10 }); + + VTKM_CONT void SetCenter(const vtkm::Vec& center) { this->Center = center; } + + VTKM_CONT void SetSpacing(const vtkm::Vec& spacing) { this->Spacing = spacing; } + + VTKM_CONT void SetFrequency(const vtkm::Vec& frequency) + { + this->Frequency = frequency; + } + + VTKM_CONT void SetMagnitude(const vtkm::Vec& magnitude) + { + this->Magnitude = magnitude; + } + + VTKM_CONT void SetMinimumExtent(const vtkm::Id3& minExtent) { this->MinimumExtent = minExtent; } + + VTKM_CONT void SetMaximumExtent(const vtkm::Id3& maxExtent) { this->MaximumExtent = maxExtent; } + + VTKM_CONT void SetExtent(const vtkm::Id3& minExtent, const vtkm::Id3& maxExtent) + { + this->MinimumExtent = minExtent; + this->MaximumExtent = maxExtent; + } + + VTKM_CONT void SetMaximumValue(const vtkm::FloatDefault& maxVal) { this->MaximumValue = maxVal; } + + VTKM_CONT void SetStandardDeviation(const vtkm::FloatDefault& stdev) + { + this->StandardDeviation = stdev; + } + + vtkm::cont::DataSet Execute() const; + +private: + vtkm::cont::Field GeneratePointField(const std::string& name) const; + + vtkm::Vec3f Center; + vtkm::Vec3f Spacing; + vtkm::Vec3f Frequency; + vtkm::Vec3f Magnitude; + vtkm::Id3 MinimumExtent; + vtkm::Id3 MaximumExtent; + vtkm::FloatDefault MaximumValue; + vtkm::FloatDefault StandardDeviation; +}; +} //namespace source +} //namespace vtkm + +#endif //vtk_m_source_Wavelet_h diff --git a/vtkm/source/testing/UnitTestWaveletSource.cxx b/vtkm/source/testing/UnitTestWaveletSource.cxx index bb5bd2410..20081e5cc 100644 --- a/vtkm/source/testing/UnitTestWaveletSource.cxx +++ b/vtkm/source/testing/UnitTestWaveletSource.cxx @@ -13,7 +13,7 @@ #include #include -void WaveletGeneratorTest() +void WaveletSourceTest() { vtkm::cont::Timer timer; timer.Start(); @@ -27,7 +27,7 @@ void WaveletGeneratorTest() std::cout << "Default wavelet took " << time << "s.\n"; { - auto coords = ds.GetCoordinateSystem("coords"); + auto coords = ds.GetCoordinateSystem("coordinates"); auto data = coords.GetData(); VTKM_TEST_ASSERT(test_equal(data.GetNumberOfValues(), 9261), "Incorrect number of points."); } @@ -63,7 +63,7 @@ void WaveletGeneratorTest() } } -int UnitTestWaveletGenerator(int argc, char* argv[]) +int UnitTestWaveletSource(int argc, char* argv[]) { - return vtkm::cont::testing::Testing::Run(WaveletGeneratorTest, argc, argv); + return vtkm::cont::testing::Testing::Run(WaveletSourceTest, argc, argv); } diff --git a/vtkm/worklet/testing/UnitTestOrientNormals.cxx b/vtkm/worklet/testing/UnitTestOrientNormals.cxx index 0ab38a96c..ce09be028 100644 --- a/vtkm/worklet/testing/UnitTestOrientNormals.cxx +++ b/vtkm/worklet/testing/UnitTestOrientNormals.cxx @@ -37,7 +37,7 @@ #include #include -#include +#include #include @@ -59,12 +59,10 @@ struct TestPolicy : public vtkm::filter::PolicyBase VTKM_CONT vtkm::cont::DataSet CreateDataSet(bool pointNormals, bool cellNormals) { - vtkm::worklet::WaveletGenerator wavelet; - wavelet.SetMinimumExtent({ -25 }); - wavelet.SetMaximumExtent({ 25 }); + vtkm::source::Wavelet wavelet({ -25 }, { 25 }); wavelet.SetFrequency({ 20, 15, 25 }); wavelet.SetMagnitude({ 5 }); - auto dataSet = wavelet.GenerateDataSet(); + auto dataSet = wavelet.Execute(); // Cut a contour vtkm::filter::Contour contour;