mirror of
https://gitlab.kitware.com/vtk/vtk-m
synced 2024-09-16 17:22:55 +00:00
ExtractPoints worklet for ImplicitFunction
Add Box to ImplicitFunction
This commit is contained in:
Patricia Kroll Fasel - 090207
parent
be73dec0cb
commit
8598711c09
@ -21,6 +21,8 @@
|
||||
#define vtk_m_ImplicitFunctions_h
|
||||
|
||||
#include <vtkm/Types.h>
|
||||
#include <vtkm/Math.h>
|
||||
#include <vtkm/VectorAnalysis.h>
|
||||
#include <iostream>
|
||||
|
||||
namespace vtkm {
|
||||
@ -153,6 +155,228 @@ private:
|
||||
vtkm::Vec<FloatDefault, 3> Center;
|
||||
};
|
||||
|
||||
/// \brief Implicit function for a box
|
||||
class Box
|
||||
{
|
||||
public:
|
||||
VTKM_CONT
|
||||
Box() : MinPoint(vtkm::Vec<FloatDefault,3>(0.0, 0.0, 0.0)),
|
||||
MaxPoint(vtkm::Vec<FloatDefault,3>(1.0, 1.0, 1.0))
|
||||
{ }
|
||||
|
||||
VTKM_CONT
|
||||
Box(vtkm::Vec<FloatDefault, 3> minPoint, vtkm::Vec<FloatDefault, 3> maxPoint)
|
||||
: MinPoint(minPoint), MaxPoint(maxPoint)
|
||||
{ }
|
||||
|
||||
VTKM_CONT
|
||||
Box(FloatDefault xmin, FloatDefault xmax,
|
||||
FloatDefault ymin, FloatDefault ymax,
|
||||
FloatDefault zmin, FloatDefault zmax)
|
||||
{
|
||||
MinPoint[0] = xmin; MaxPoint[0] = xmax;
|
||||
MinPoint[1] = ymin; MaxPoint[1] = ymax;
|
||||
MinPoint[2] = zmin; MaxPoint[2] = zmax;
|
||||
}
|
||||
|
||||
VTKM_EXEC_CONT
|
||||
const vtkm::Vec<FloatDefault, 3>& GetMinPoint() const
|
||||
{
|
||||
return this->MinPoint;
|
||||
}
|
||||
|
||||
VTKM_EXEC_CONT
|
||||
const vtkm::Vec<FloatDefault, 3>& GetMaxPoint() const
|
||||
{
|
||||
return this->MaxPoint;
|
||||
}
|
||||
|
||||
VTKM_EXEC_CONT
|
||||
FloatDefault Value(const vtkm::Vec<FloatDefault, 3> &x) const
|
||||
{
|
||||
vtkm::Vec<vtkm::IdComponent,3> inside(1, 1, 1);
|
||||
vtkm::Vec<FloatDefault,3> dist(0.0, 0.0, 0.0);
|
||||
FloatDefault insideDistance = vtkm::NegativeInfinity32();
|
||||
|
||||
for (vtkm::IdComponent d = 0; d < 3; d++)
|
||||
{
|
||||
if (this->MinPoint[d] == this->MaxPoint[d])
|
||||
{
|
||||
dist[d] = vtkm::Abs(x[d] - MinPoint[d]);
|
||||
if (dist[d] > 0.0)
|
||||
inside[d] = 0;
|
||||
}
|
||||
else
|
||||
{
|
||||
// Calculate the distance of point to box boundary
|
||||
if (x[d] < this->MinPoint[d])
|
||||
{
|
||||
// Point less than bounding box minimum (positive dist)
|
||||
inside[d] = 0;
|
||||
dist[d] = this->MinPoint[d] - x[d];
|
||||
}
|
||||
else if (x[d] > this->MaxPoint[d])
|
||||
{
|
||||
// Point greater than bounding box maximum (positive dist)
|
||||
inside[d] = 0;
|
||||
dist[d] = x[d] - this->MaxPoint[d];
|
||||
}
|
||||
else if (x[d] <= ((this->MaxPoint[d] - this->MinPoint[d]) / 2.0))
|
||||
{
|
||||
// Point inside box closer to minimum (negative dist)
|
||||
dist[d] = this->MinPoint[d] - x[d];
|
||||
if (dist[d] > insideDistance)
|
||||
insideDistance = dist[d];
|
||||
}
|
||||
else
|
||||
{
|
||||
// Point inside box closer to maximum (negative dist)
|
||||
dist[d] = x[d] - this->MaxPoint[d];
|
||||
if (dist[d] > insideDistance)
|
||||
insideDistance = dist[d];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (inside[0] && inside[1] && inside[2])
|
||||
{
|
||||
return(insideDistance);
|
||||
}
|
||||
else
|
||||
{
|
||||
FloatDefault distance = 0.0;
|
||||
for (vtkm::IdComponent d = 0; d < 3; d++)
|
||||
{
|
||||
if (dist[d] > 0.0)
|
||||
distance += dist[d] * dist[d];
|
||||
}
|
||||
return vtkm::Sqrt(distance);
|
||||
}
|
||||
}
|
||||
|
||||
VTKM_EXEC_CONT
|
||||
FloatDefault Value(FloatDefault x, FloatDefault y, FloatDefault z) const
|
||||
{
|
||||
return this->Value(vtkm::Vec<vtkm::FloatDefault,3>(x, y, z));
|
||||
}
|
||||
|
||||
VTKM_EXEC_CONT
|
||||
vtkm::Vec<FloatDefault, 3> Gradient(const vtkm::Vec<FloatDefault, 3> &x) const
|
||||
{
|
||||
vtkm::IdComponent minAxis = 0;
|
||||
FloatDefault dist = 0.0;
|
||||
FloatDefault minDist = vtkm::Infinity32();
|
||||
vtkm::Vec<vtkm::IdComponent,3> location;
|
||||
vtkm::Vec<FloatDefault,3> normal;
|
||||
vtkm::Vec<FloatDefault,3> inside(0.0, 0.0, 0.0);
|
||||
vtkm::Vec<FloatDefault,3> outside(0.0, 0.0, 0.0);
|
||||
vtkm::Vec<FloatDefault,3> center((this->MaxPoint[0] - this->MinPoint[0]) / 2.0,
|
||||
(this->MaxPoint[1] - this->MinPoint[1]) / 2.0,
|
||||
(this->MaxPoint[2] - this->MinPoint[2]) / 2.0);
|
||||
|
||||
// Compute the location of the point with respect to the box
|
||||
// Point will lie in one of 27 separate regions around or within the box
|
||||
// Gradient vector is computed differently in each of the regions.
|
||||
for (vtkm::IdComponent d = 0; d < 3; d++)
|
||||
{
|
||||
if (x[d] < this->MinPoint[d])
|
||||
{
|
||||
// Outside the box low end
|
||||
location[d] = 0;
|
||||
outside[d] = -1.0;
|
||||
}
|
||||
else if (x[d] > this->MaxPoint[d])
|
||||
{
|
||||
// Outside the box high end
|
||||
location[d] = 2;
|
||||
outside[d] = 1.0;
|
||||
}
|
||||
else
|
||||
{
|
||||
if (x[d] <= center[d])
|
||||
{
|
||||
// Inside the box low end
|
||||
location[d] = 1;
|
||||
inside[d] = -1.0;
|
||||
dist = x[d] - this->MinPoint[d];
|
||||
}
|
||||
else
|
||||
{
|
||||
// Inside the box high end
|
||||
location[d] = 1;
|
||||
inside[d] = 1.0;
|
||||
dist = this->MaxPoint[d] - x[d];
|
||||
}
|
||||
if (dist < minDist) // dist is negative
|
||||
{
|
||||
minDist = dist;
|
||||
minAxis = d;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
int indx = location[0] + 3*location[1] + 9*location[2];
|
||||
switch (indx)
|
||||
{
|
||||
// verts - gradient points away from center point
|
||||
case 0: case 2: case 6: case 8: case 18: case 20: case 24: case 26:
|
||||
for (vtkm::IdComponent d = 0; d < 3; d++)
|
||||
{
|
||||
normal[d] = x[d] - center[d];
|
||||
}
|
||||
vtkm::Normalize(normal);
|
||||
break;
|
||||
|
||||
// edges - gradient points out from axis of cube
|
||||
case 1: case 3: case 5: case 7:
|
||||
case 9: case 11: case 15: case 17:
|
||||
case 19: case 21: case 23: case 25:
|
||||
for (vtkm::IdComponent d = 0; d < 3; d++)
|
||||
{
|
||||
if (outside[d] != 0.0)
|
||||
{
|
||||
normal[d] = x[d] - center[d];
|
||||
}
|
||||
else
|
||||
{
|
||||
normal[d] = 0.0;
|
||||
}
|
||||
}
|
||||
vtkm::Normalize(normal);
|
||||
break;
|
||||
|
||||
// faces - gradient points perpendicular to face
|
||||
case 4: case 10: case 12: case 14: case 16: case 22:
|
||||
for (vtkm::IdComponent d = 0; d < 3; d++)
|
||||
{
|
||||
normal[d] = outside[d];
|
||||
}
|
||||
break;
|
||||
|
||||
// interior - gradient is perpendicular to closest face
|
||||
case 13:
|
||||
normal[0] = normal[1] = normal[2] = 0.0;
|
||||
normal[minAxis] = inside[minAxis];
|
||||
break;
|
||||
default:
|
||||
assert("check: impossible case." && 0); // reaching this line is a bug.
|
||||
break;
|
||||
}
|
||||
return normal;
|
||||
}
|
||||
|
||||
VTKM_EXEC_CONT
|
||||
vtkm::Vec<FloatDefault, 3> Gradient(FloatDefault x, FloatDefault y, FloatDefault z)
|
||||
const
|
||||
{
|
||||
return this->Gradient(vtkm::Vec<FloatDefault, 3>(x, y, z));
|
||||
}
|
||||
|
||||
private:
|
||||
vtkm::Vec<FloatDefault, 3> MinPoint;
|
||||
vtkm::Vec<FloatDefault, 3> MaxPoint;
|
||||
};
|
||||
|
||||
/// \brief A function object that evaluates the contained implicit function
|
||||
template <typename ImplicitFunction>
|
||||
class ImplicitFunctionValue
|
||||
|
||||
@ -29,6 +29,7 @@ set(headers
|
||||
DispatcherReduceByKey.h
|
||||
DispatcherStreamingMapField.h
|
||||
ExternalFaces.h
|
||||
ExtractPoints.h
|
||||
FieldHistogram.h
|
||||
FieldStatistics.h
|
||||
Gradient.h
|
||||
|
||||
108
vtkm/worklet/ExtractPoints.h
Normal file
108
vtkm/worklet/ExtractPoints.h
Normal file
@ -0,0 +1,108 @@
|
||||
//============================================================================
|
||||
// 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.
|
||||
//
|
||||
// Copyright 2014 Sandia Corporation.
|
||||
// Copyright 2014 UT-Battelle, LLC.
|
||||
// Copyright 2014 Los Alamos National Security.
|
||||
//
|
||||
// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
|
||||
// the U.S. Government retains certain rights in this software.
|
||||
//
|
||||
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
|
||||
// Laboratory (LANL), the U.S. Government retains certain rights in
|
||||
// this software.
|
||||
//============================================================================
|
||||
#ifndef vtkm_m_worklet_ExtractPoints_h
|
||||
#define vtkm_m_worklet_ExtractPoints_h
|
||||
|
||||
#include <vtkm/worklet/ScatterCounting.h>
|
||||
#include <vtkm/worklet/DispatcherMapTopology.h>
|
||||
#include <vtkm/worklet/WorkletMapTopology.h>
|
||||
#include <vtkm/worklet/WorkletMapField.h>
|
||||
|
||||
#include <vtkm/cont/DataSet.h>
|
||||
#include <vtkm/cont/ArrayHandle.h>
|
||||
#include <vtkm/cont/CoordinateSystem.h>
|
||||
#include <vtkm/cont/DeviceAdapterAlgorithm.h>
|
||||
#include <vtkm/ImplicitFunctions.h>
|
||||
|
||||
namespace vtkm {
|
||||
namespace worklet {
|
||||
|
||||
class ExtractPoints : public vtkm::worklet::WorkletMapPointToCell
|
||||
{
|
||||
public:
|
||||
ExtractPoints() {}
|
||||
|
||||
template<typename ImplicitFunction>
|
||||
class ExtractPointsWithImplicitFunction : public vtkm::worklet::WorkletMapField
|
||||
{
|
||||
public:
|
||||
typedef void ControlSignature(FieldIn<Vec3> coordinates,
|
||||
FieldOut<IdComponentType> mask);
|
||||
typedef _2 ExecutionSignature(_1);
|
||||
|
||||
VTKM_CONT
|
||||
ExtractPointsWithImplicitFunction(const ImplicitFunction &function) :
|
||||
Function(function) {}
|
||||
|
||||
VTKM_CONT
|
||||
vtkm::IdComponent operator()(const vtkm::Vec<vtkm::Float64,3> &coordinate) const
|
||||
{
|
||||
vtkm::Float64 value = this->Function.Value(coordinate);
|
||||
vtkm::Float64 gradient = this->Function.Value(coordinate);
|
||||
vtkm::IdComponent mask = 0;
|
||||
if (value <= 0)
|
||||
mask = 1;
|
||||
std::cout << "Coord " << coordinate[0] << " , " << coordinate[1] << " , " << coordinate[2] << " value " << value << " mask " << mask << " GRADIENT " << gradient << std::endl;
|
||||
return mask;
|
||||
}
|
||||
|
||||
private:
|
||||
ImplicitFunction Function;
|
||||
};
|
||||
|
||||
template <typename CellSetType,
|
||||
typename ImplicitFunction,
|
||||
typename DeviceAdapter>
|
||||
vtkm::cont::CellSetSingleType<> Run(
|
||||
const CellSetType &cellSet,
|
||||
const ImplicitFunction &implicitFunction,
|
||||
const vtkm::cont::CoordinateSystem &coordinates,
|
||||
DeviceAdapter device)
|
||||
{
|
||||
typedef typename vtkm::cont::DeviceAdapterAlgorithm<DeviceAdapter> DeviceAlgorithm;
|
||||
vtkm::cont::ArrayHandle<vtkm::IdComponent> maskArray;
|
||||
vtkm::Id numberOfInputPoints = cellSet.GetNumberOfPoints();
|
||||
DeviceAlgorithm::Copy(vtkm::cont::ArrayHandleConstant<vtkm::IdComponent>(0, numberOfInputPoints),
|
||||
maskArray);
|
||||
|
||||
// Worklet output will be a boolean passFlag array
|
||||
typedef ExtractPointsWithImplicitFunction<ImplicitFunction> ExtractPointsWorklet;
|
||||
ExtractPointsWorklet worklet(implicitFunction);
|
||||
DispatcherMapField<ExtractPointsWorklet, DeviceAdapter> dispatcher(worklet);
|
||||
dispatcher.Invoke(coordinates, maskArray);
|
||||
|
||||
vtkm::worklet::ScatterCounting PointScatter(maskArray, DeviceAdapter(), true);
|
||||
vtkm::cont::ArrayHandle<vtkm::Id> pointIds = PointScatter.GetOutputToInputMap();
|
||||
|
||||
// Make CellSetSingleType with VERTEX at each point id
|
||||
vtkm::cont::CellSetSingleType< > outCellSet("cells");
|
||||
outCellSet.Fill(numberOfInputPoints,
|
||||
vtkm::CellShapeTagVertex::Id,
|
||||
1,
|
||||
pointIds);
|
||||
|
||||
return outCellSet;
|
||||
}
|
||||
};
|
||||
|
||||
}
|
||||
} // namespace vtkm::worklet
|
||||
|
||||
#endif // vtkm_m_worklet_ExtractPoints_h
|
||||
@ -45,7 +45,6 @@ public:
|
||||
public:
|
||||
typedef void ControlSignature(FieldIn<> scalars,
|
||||
FieldOut<IdComponentType> mask);
|
||||
|
||||
typedef _2 ExecutionSignature(_1);
|
||||
|
||||
VTKM_CONT
|
||||
|
||||
@ -25,6 +25,7 @@ set(unit_tests
|
||||
UnitTestClipping.cxx
|
||||
UnitTestContourTreeUniform.cxx
|
||||
UnitTestExternalFaces.cxx
|
||||
UnitTestExtractPoints.cxx
|
||||
UnitTestFieldHistogram.cxx
|
||||
UnitTestFieldStatistics.cxx
|
||||
UnitTestKeys.cxx
|
||||
|
||||
133
vtkm/worklet/testing/UnitTestExtractPoints.cxx
Normal file
133
vtkm/worklet/testing/UnitTestExtractPoints.cxx
Normal file
@ -0,0 +1,133 @@
|
||||
//============================================================================
|
||||
// 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.
|
||||
//
|
||||
// Copyright 2014 Sandia Corporation.
|
||||
// Copyright 2014 UT-Battelle, LLC.
|
||||
// Copyright 2014 Los Alamos National Security.
|
||||
//
|
||||
// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
|
||||
// the U.S. Government retains certain rights in this software.
|
||||
//
|
||||
// Under the terms of Contract DE-AC52-06NA25396 with Los Alamos National
|
||||
// Laboratory (LANL), the U.S. Government retains certain rights in
|
||||
// this software.
|
||||
//============================================================================
|
||||
|
||||
#include <vtkm/worklet/ExtractPoints.h>
|
||||
|
||||
#include <vtkm/cont/testing/MakeTestDataSet.h>
|
||||
#include <vtkm/cont/testing/Testing.h>
|
||||
|
||||
#include <vtkm/cont/ArrayPortalToIterators.h>
|
||||
#include <vtkm/cont/CellSet.h>
|
||||
|
||||
#include <algorithm>
|
||||
#include <iostream>
|
||||
#include <vector>
|
||||
|
||||
using vtkm::cont::testing::MakeTestDataSet;
|
||||
|
||||
template <typename DeviceAdapter>
|
||||
class TestingExtractPoints
|
||||
{
|
||||
public:
|
||||
void TestExtractPointsWithSphere() const
|
||||
{
|
||||
std::cout << "Testing extract points with implicit function (sphere):" << std::endl;
|
||||
|
||||
typedef vtkm::cont::CellSetStructured<3> CellSetType;
|
||||
typedef vtkm::cont::CellSetSingleType<> OutCellSetType;
|
||||
|
||||
// Input data set created
|
||||
vtkm::cont::DataSet dataset = MakeTestDataSet().Make3DUniformDataSet1();
|
||||
CellSetType cellset;
|
||||
dataset.GetCellSet(0).CopyTo(cellset);
|
||||
|
||||
// Implicit function
|
||||
vtkm::Vec<vtkm::FloatDefault, 3> center(2, 2, 2);
|
||||
vtkm::FloatDefault radius(1.8);
|
||||
vtkm::Sphere sphere(center, radius);
|
||||
|
||||
// Output dataset contains input coordinate system and point data
|
||||
vtkm::cont::DataSet outDataSet;
|
||||
outDataSet.AddCoordinateSystem(dataset.GetCoordinateSystem(0));
|
||||
for (vtkm::Id indx = 0; indx < dataset.GetNumberOfFields(); indx++)
|
||||
{
|
||||
vtkm::cont::Field field = dataset.GetField(indx);
|
||||
if (field.GetAssociation() == vtkm::cont::Field::ASSOC_POINTS)
|
||||
{
|
||||
outDataSet.AddField(field);
|
||||
}
|
||||
}
|
||||
|
||||
// Output data set with cell set containing extracted points
|
||||
vtkm::worklet::ExtractPoints extractPoints;
|
||||
OutCellSetType outCellSet;
|
||||
vtkm::cont::CellSetSingleType<> outputCellSet =
|
||||
outCellSet = extractPoints.Run(dataset.GetCellSet(0),
|
||||
sphere,
|
||||
dataset.GetCoordinateSystem("coords"),
|
||||
DeviceAdapter());
|
||||
outDataSet.AddCellSet(outCellSet);
|
||||
|
||||
VTKM_TEST_ASSERT(test_equal(outCellSet.GetNumberOfCells(), 27), "Wrong result for ExtractPoints");
|
||||
}
|
||||
|
||||
void TestExtractPointsWithBox() const
|
||||
{
|
||||
std::cout << "Testing extract points with implicit function (box):" << std::endl;
|
||||
|
||||
typedef vtkm::cont::CellSetStructured<3> CellSetType;
|
||||
typedef vtkm::cont::CellSetSingleType<> OutCellSetType;
|
||||
|
||||
// Input data set created
|
||||
vtkm::cont::DataSet dataset = MakeTestDataSet().Make3DUniformDataSet1();
|
||||
CellSetType cellset;
|
||||
dataset.GetCellSet(0).CopyTo(cellset);
|
||||
|
||||
// Implicit function
|
||||
vtkm::Vec<vtkm::FloatDefault, 3> minPoint(1.0, 1.0, 1.0);
|
||||
vtkm::Vec<vtkm::FloatDefault, 3> maxPoint(3.0, 3.0, 3.0);
|
||||
vtkm::Box box(minPoint, maxPoint);
|
||||
|
||||
// Output dataset contains input coordinate system and point data
|
||||
vtkm::cont::DataSet outDataSet;
|
||||
outDataSet.AddCoordinateSystem(dataset.GetCoordinateSystem(0));
|
||||
for (vtkm::Id indx = 0; indx < dataset.GetNumberOfFields(); indx++)
|
||||
{
|
||||
vtkm::cont::Field field = dataset.GetField(indx);
|
||||
if (field.GetAssociation() == vtkm::cont::Field::ASSOC_POINTS)
|
||||
{
|
||||
outDataSet.AddField(field);
|
||||
}
|
||||
}
|
||||
|
||||
// Output data set with cell set containing extracted points
|
||||
vtkm::worklet::ExtractPoints extractPoints;
|
||||
OutCellSetType outCellSet;
|
||||
vtkm::cont::CellSetSingleType<> outputCellSet =
|
||||
outCellSet = extractPoints.Run(dataset.GetCellSet(0),
|
||||
box,
|
||||
dataset.GetCoordinateSystem("coords"),
|
||||
DeviceAdapter());
|
||||
outDataSet.AddCellSet(outCellSet);
|
||||
|
||||
VTKM_TEST_ASSERT(test_equal(outCellSet.GetNumberOfCells(), 27), "Wrong result for ExtractPoints");
|
||||
}
|
||||
void operator()() const
|
||||
{
|
||||
this->TestExtractPointsWithSphere();
|
||||
this->TestExtractPointsWithBox();
|
||||
}
|
||||
};
|
||||
|
||||
int UnitTestExtractPoints(int, char *[])
|
||||
{
|
||||
return vtkm::cont::testing::Testing::Run(
|
||||
TestingExtractPoints<VTKM_DEFAULT_DEVICE_ADAPTER_TAG>());
|
||||
}
|
||||
Loading…
Reference in New Issue
Block a user