Merge topic 'flying_edges_support_multi_contour'

243425421 Flying Edges now works with multiple contours

Acked-by: Kitware Robot <kwrobot@kitware.com>
Acked-by: Kenneth Moreland <kmorel@sandia.gov>
Merge-request: !2039
This commit is contained in:
Robert Maynard 2020-04-15 17:46:16 +00:00 committed by Kitware Robot
commit fdffaac2d2
3 changed files with 71 additions and 33 deletions

@ -128,6 +128,8 @@ vtkm::cont::CellSetSingleType<> execute(
vtkm::cont::ArrayHandle<vtkm::Id> triangle_topology; vtkm::cont::ArrayHandle<vtkm::Id> triangle_topology;
for (std::size_t i = 0; i < isovalues.size(); ++i) for (std::size_t i = 0; i < isovalues.size(); ++i)
{ {
auto multiContourCellOffset = sharedState.CellIdMap.GetNumberOfValues();
auto multiContourPointOffset = sharedState.InterpolationWeights.GetNumberOfValues();
ValueType isoval = isovalues[i]; ValueType isoval = isovalues[i];
//---------------------------------------------------------------------------- //----------------------------------------------------------------------------
@ -161,16 +163,16 @@ vtkm::cont::CellSetSingleType<> execute(
detail::extend_by(triangle_topology, 3 * sumTris); detail::extend_by(triangle_topology, 3 * sumTris);
detail::extend_by(sharedState.CellIdMap, sumTris); detail::extend_by(sharedState.CellIdMap, sumTris);
//By using the current length of points array as the starting value we
//handle multiple contours, by propagating the correct write offset auto newPointSize =
auto newPointSize = vtkm::cont::Algorithm::ScanExclusive( vtkm::cont::Algorithm::ScanExclusive(metaDataLinearSums, metaDataLinearSums);
metaDataLinearSums, metaDataLinearSums, vtkm::Sum{}, points.GetNumberOfValues());
detail::extend_by(sharedState.InterpolationEdgeIds, newPointSize); detail::extend_by(sharedState.InterpolationEdgeIds, newPointSize);
detail::extend_by(sharedState.InterpolationWeights, newPointSize); detail::extend_by(sharedState.InterpolationWeights, newPointSize);
//---------------------------------------------------------------------------- //----------------------------------------------------------------------------
// PASS 4: Process voxel rows and generate topology, and interpolation state // PASS 4: Process voxel rows and generate topology, and interpolation state
ComputePass4<ValueType, AxisToSum> worklet4(isoval, pdims); ComputePass4<ValueType, AxisToSum> worklet4(
isoval, pdims, multiContourCellOffset, multiContourPointOffset);
invoke(worklet4, invoke(worklet4,
metaDataMesh2D, metaDataMesh2D,
metaDataSums, metaDataSums,

@ -55,7 +55,7 @@ VTKM_EXEC inline void generate_tris(vtkm::Id inputCellId,
vtkm::UInt8 edgeCase, vtkm::UInt8 edgeCase,
vtkm::UInt8 numTris, vtkm::UInt8 numTris,
vtkm::Id* edgeIds, vtkm::Id* edgeIds,
vtkm::Int32& triId, vtkm::Id& triId,
const WholeConnField& conn, const WholeConnField& conn,
const WholeCellIdField& cellIds) const WholeCellIdField& cellIds)
{ {
@ -84,22 +84,23 @@ VTKM_EXEC inline void generate_tris(vtkm::Id inputCellId,
//---------------------------------------------------------------------------- //----------------------------------------------------------------------------
template <typename AxisToSum, typename FieldInPointId3> template <typename AxisToSum, typename FieldInPointId3>
VTKM_EXEC inline void init_voxelIds(AxisToSum, VTKM_EXEC inline void init_voxelIds(AxisToSum,
vtkm::Id writeOffset,
vtkm::UInt8 edgeCase, vtkm::UInt8 edgeCase,
const FieldInPointId3& axis_sums, const FieldInPointId3& axis_sums,
vtkm::Id* edgeIds) vtkm::Id* edgeIds)
{ {
auto* edgeUses = data::GetEdgeUses(edgeCase); auto* edgeUses = data::GetEdgeUses(edgeCase);
edgeIds[0] = axis_sums[0][AxisToSum::xindex]; // x-edges edgeIds[0] = writeOffset + axis_sums[0][AxisToSum::xindex]; // x-edges
edgeIds[1] = axis_sums[1][AxisToSum::xindex]; edgeIds[1] = writeOffset + axis_sums[1][AxisToSum::xindex];
edgeIds[2] = axis_sums[3][AxisToSum::xindex]; edgeIds[2] = writeOffset + axis_sums[3][AxisToSum::xindex];
edgeIds[3] = axis_sums[2][AxisToSum::xindex]; edgeIds[3] = writeOffset + axis_sums[2][AxisToSum::xindex];
edgeIds[4] = axis_sums[0][AxisToSum::yindex]; // y-edges edgeIds[4] = writeOffset + axis_sums[0][AxisToSum::yindex]; // y-edges
edgeIds[5] = edgeIds[4] + edgeUses[4]; edgeIds[5] = edgeIds[4] + edgeUses[4];
edgeIds[6] = axis_sums[3][AxisToSum::yindex]; edgeIds[6] = writeOffset + axis_sums[3][AxisToSum::yindex];
edgeIds[7] = edgeIds[6] + edgeUses[6]; edgeIds[7] = edgeIds[6] + edgeUses[6];
edgeIds[8] = axis_sums[0][AxisToSum::zindex]; // z-edges edgeIds[8] = writeOffset + axis_sums[0][AxisToSum::zindex]; // z-edges
edgeIds[9] = edgeIds[8] + edgeUses[8]; edgeIds[9] = edgeIds[8] + edgeUses[8];
edgeIds[10] = axis_sums[1][AxisToSum::zindex]; edgeIds[10] = writeOffset + axis_sums[1][AxisToSum::zindex];
edgeIds[11] = edgeIds[10] + edgeUses[10]; edgeIds[11] = edgeIds[10] + edgeUses[10];
} }
@ -121,19 +122,24 @@ VTKM_EXEC inline void advance_voxelIds(vtkm::UInt8 const* const edgeUses, vtkm::
edgeIds[11] = edgeIds[10] + edgeUses[11]; edgeIds[11] = edgeIds[10] + edgeUses[11];
} }
template <typename T, typename AxisToSum> template <typename T, typename AxisToSum>
struct ComputePass4 : public vtkm::worklet::WorkletVisitCellsWithPoints struct ComputePass4 : public vtkm::worklet::WorkletVisitCellsWithPoints
{ {
vtkm::Id3 PointDims; vtkm::Id3 PointDims;
T IsoValue; T IsoValue;
vtkm::Id CellWriteOffset;
vtkm::Id PointWriteOffset;
ComputePass4() {} ComputePass4() {}
ComputePass4(T value, const vtkm::Id3& pdims) ComputePass4(T value,
const vtkm::Id3& pdims,
vtkm::Id multiContourCellOffset,
vtkm::Id multiContourPointOffset)
: PointDims(pdims) : PointDims(pdims)
, IsoValue(value) , IsoValue(value)
, CellWriteOffset(multiContourCellOffset)
, PointWriteOffset(multiContourPointOffset)
{ {
} }
@ -176,12 +182,13 @@ struct ComputePass4 : public vtkm::worklet::WorkletVisitCellsWithPoints
{ {
//This works as cellTriCount was computed with ScanExtended //This works as cellTriCount was computed with ScanExtended
//and therefore has one more entry than the number of cells //and therefore has one more entry than the number of cells
vtkm::Int32 cell_tri_offset = cellTriCount.Get(oidx); vtkm::Id cell_tri_offset = cellTriCount.Get(oidx);
vtkm::Int32 next_tri_offset = cellTriCount.Get(oidx + 1); vtkm::Id next_tri_offset = cellTriCount.Get(oidx + 1);
if (cell_tri_offset == next_tri_offset) if (cell_tri_offset == next_tri_offset)
{ //we produce nothing { //we produce nothing
return; return;
} }
cell_tri_offset += this->CellWriteOffset;
// find adjusted trim values. // find adjusted trim values.
vtkm::Id left = vtkm::Min(axis_mins[0], axis_mins[1]); vtkm::Id left = vtkm::Min(axis_mins[0], axis_mins[1]);
@ -242,7 +249,7 @@ struct ComputePass4 : public vtkm::worklet::WorkletVisitCellsWithPoints
const vtkm::Id3 increments = compute_incs3d(pdims); const vtkm::Id3 increments = compute_incs3d(pdims);
vtkm::Id edgeIds[12]; vtkm::Id edgeIds[12];
auto edgeCase = getEdgeCase(edges, startPos, (axis_inc * left)); auto edgeCase = getEdgeCase(edges, startPos, (axis_inc * left));
init_voxelIds(AxisToSum{}, edgeCase, axis_sums, edgeIds); init_voxelIds(AxisToSum{}, this->PointWriteOffset, edgeCase, axis_sums, edgeIds);
for (vtkm::Id i = left; i < right; ++i) // run along the trimmed voxels for (vtkm::Id i = left; i < right; ++i) // run along the trimmed voxels
{ {
ijk[AxisToSum::xindex] = i; ijk[AxisToSum::xindex] = i;

@ -196,25 +196,26 @@ void TestContourUniformGrid()
vtkm::cont::ArrayHandle<vtkm::FloatDefault> cellFieldArray; vtkm::cont::ArrayHandle<vtkm::FloatDefault> cellFieldArray;
dataSet.GetField("cellvar").GetData().CopyTo(cellFieldArray); dataSet.GetField("cellvar").GetData().CopyTo(cellFieldArray);
vtkm::worklet::Contour isosurfaceFilter; vtkm::worklet::Contour contour;
isosurfaceFilter.SetMergeDuplicatePoints(false); contour.SetMergeDuplicatePoints(false);
std::vector<vtkm::Float32> contourValue{ 0.5f }; std::vector<vtkm::Float32> contourValue{ 0.5f, 0.5f };
const vtkm::Id numContours = static_cast<vtkm::Id>(contourValue.size());
vtkm::cont::ArrayHandle<vtkm::Vec3f_32> verticesArray; vtkm::cont::ArrayHandle<vtkm::Vec3f_32> verticesArray;
vtkm::cont::ArrayHandle<vtkm::Vec3f_32> normalsArray; vtkm::cont::ArrayHandle<vtkm::Vec3f_32> normalsArray;
vtkm::cont::ArrayHandle<vtkm::Float32> scalarsArray; vtkm::cont::ArrayHandle<vtkm::Float32> scalarsArray;
auto result = isosurfaceFilter.Run(contourValue, auto result = contour.Run(contourValue,
cellSet, cellSet,
dataSet.GetCoordinateSystem(), dataSet.GetCoordinateSystem(),
pointFieldArray, pointFieldArray,
verticesArray, verticesArray,
normalsArray); normalsArray);
scalarsArray = isosurfaceFilter.ProcessPointField(pointFieldArray); scalarsArray = contour.ProcessPointField(pointFieldArray);
vtkm::cont::ArrayHandle<vtkm::FloatDefault> cellFieldArrayOut; vtkm::cont::ArrayHandle<vtkm::FloatDefault> cellFieldArrayOut;
cellFieldArrayOut = isosurfaceFilter.ProcessCellField(cellFieldArray); cellFieldArrayOut = contour.ProcessCellField(cellFieldArray);
std::cout << "vertices: "; std::cout << "vertices: ";
vtkm::cont::printSummary_ArrayHandle(verticesArray, std::cout); vtkm::cont::printSummary_ArrayHandle(verticesArray, std::cout);
@ -231,9 +232,37 @@ void TestContourUniformGrid()
VTKM_TEST_ASSERT(result.GetNumberOfCells() == cellFieldArrayOut.GetNumberOfValues()); VTKM_TEST_ASSERT(result.GetNumberOfCells() == cellFieldArrayOut.GetNumberOfValues());
VTKM_TEST_ASSERT(result.GetNumberOfCells() == 160); VTKM_TEST_ASSERT(result.GetNumberOfCells() == (160 * numContours));
VTKM_TEST_ASSERT(verticesArray.GetNumberOfValues() == 72); VTKM_TEST_ASSERT(verticesArray.GetNumberOfValues() == (72 * numContours));
// Verify that multiple contours of the same iso value are identical
{
auto normal_portal = normalsArray.ReadPortal();
for (vtkm::Id i = 0; i < 72; ++i)
{
for (vtkm::Id j = 1; j < numContours; ++j)
{
vtkm::Id jIndex = i + (72 * j);
VTKM_TEST_ASSERT(test_equal(normal_portal.Get(i), normal_portal.Get(jIndex)),
"multi contour failed");
}
}
auto outCellPortal =
result.GetConnectivityArray(vtkm::TopologyElementTagCell(), vtkm::TopologyElementTagPoint())
.ReadPortal();
for (vtkm::Id i = 0; i < 480; ++i)
{ //(3*160) as we are iterating triangle soup so the length is numOfCells*3
for (vtkm::Id j = 1; j < numContours; ++j)
{
vtkm::Id jIndex = i + (480 * j);
vtkm::Id expectedValue = (72 * j) + outCellPortal.Get(i);
VTKM_TEST_ASSERT(test_equal(outCellPortal.Get(jIndex), expectedValue),
"multi contour failed");
}
}
}
} }
void TestContourExplicit() void TestContourExplicit()