Compute volume and select top volume branches

This commit is contained in:
Mingzhe Li 2023-08-30 00:46:25 -06:00 committed by Gunther H. Weber
parent 3406ba3113
commit a77e5b3c63
24 changed files with 1957 additions and 51 deletions

@ -74,6 +74,7 @@
#include <vtkm/filter/scalar_topology/ContourTreeUniformDistributed.h>
#include <vtkm/filter/scalar_topology/DistributedBranchDecompositionFilter.h>
#include <vtkm/filter/scalar_topology/SelectTopVolumeContoursFilter.h>
#include <vtkm/filter/scalar_topology/worklet/branch_decomposition/HierarchicalVolumetricBranchDecomposer.h>
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/PrintVectors.h>
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/ProcessContourTree.h>
@ -213,8 +214,8 @@ int main(int argc, char* argv[])
if (!augmentHierarchicalTree)
{
VTKM_LOG_S(vtkm::cont::LogLevel::Warn,
"Warning: --computeVolumeBranchDecomposition only "
"allowed augmentation. Enabling --augmentHierarchicalTree option.");
"Warning: --computeVolumeBranchDecomposition requires augmentation. "
"Enabling --augmentHierarchicalTree option.");
augmentHierarchicalTree = true;
}
}
@ -265,6 +266,31 @@ int main(int argc, char* argv[])
}
}
vtkm::Id numBranches = 0;
if (parser.hasOption("--numBranches"))
{
numBranches = std::stoi(parser.getOption("--numBranches"));
if (!computeHierarchicalVolumetricBranchDecomposition)
{
VTKM_LOG_S(vtkm::cont::LogLevel::Warn,
"Warning: --numBranches requires computing branch decomposition. "
"Enabling --computeHierarchicalVolumetricBranchDecomposition option.");
computeHierarchicalVolumetricBranchDecomposition = true;
}
if (!augmentHierarchicalTree)
{
VTKM_LOG_S(vtkm::cont::LogLevel::Warn,
"Warning: --numBranches requires augmentation. "
"Enabling --augmentHierarchicalTree option.");
augmentHierarchicalTree = true;
}
}
ValueType eps = 0.00001f;
if (parser.hasOption("--eps"))
eps = std::stof(parser.getOption("--eps"));
#ifdef ENABLE_HDFIO
std::string dataset_name = "data";
if (parser.hasOption("--dataset"))
@ -336,7 +362,12 @@ int main(int argc, char* argv[])
std::cout << "--augmentHierarchicalTree Augment the hierarchical tree." << std::endl;
std::cout << "--computeVolumeBranchDecomposition Compute the volume branch decomposition. "
<< std::endl;
std::cout << " Requries --augmentHierarchicalTree to be set." << std::endl;
std::cout << " Requires --augmentHierarchicalTree to be set." << std::endl;
std::cout << "--numBranches Number of top volume branches to select." << std::endl;
std::cout << " Requires --computeVolumeBranchDecomposition." << std::endl;
std::cout
<< "--eps=<float> Floating point offset awary from the critical point. (default=0.00001)"
<< std::endl;
std::cout << "--preSplitFiles Input data is already pre-split into blocks." << std::endl;
std::cout << "--saveDot Save DOT files of the distributed contour tree " << std::endl
<< " computation (Default=False). " << std::endl;
@ -383,6 +414,8 @@ int main(int argc, char* argv[])
<< " saveOutputData=" << saveOutputData << std::endl
<< " forwardSummary=" << forwardSummary << std::endl
<< " nblocks=" << numBlocks << std::endl
<< " nbranches=" << numBranches << std::endl
<< " eps=" << eps << std::endl
#ifdef ENABLE_HDFIO
<< " dataset=" << dataset_name << " (HDF5 only)" << std::endl
<< " blocksPerDim=" << blocksPerDimIn[0] << "," << blocksPerDimIn[1] << ","
@ -642,6 +675,15 @@ int main(int argc, char* argv[])
vtkm::Float64 branchDecompTime = currTime - prevTime;
prevTime = currTime;
// Compute SelectTopVolumeContours if needed
vtkm::cont::PartitionedDataSet tp_result;
if (numBranches > 0)
{
vtkm::filter::scalar_topology::SelectTopVolumeContoursFilter tp_filter;
tp_filter.SetSavedBranches(numBranches);
tp_result = tp_filter.Execute(bd_result);
}
// Save output
if (saveOutputData)
{
@ -684,6 +726,66 @@ int main(int argc, char* argv[])
<< lowerEndGRId.Get(branch) << std::endl;
}
}
if (numBranches > 0)
{
#ifndef DEBUG_PRINT
bool print_to_files = (rank == 0);
vtkm::Id max_blocks_to_print = 1;
#else
bool print_to_files = true;
vtkm::Id max_blocks_to_print = result.GetNumberOfPartitions();
#endif
for (vtkm::Id ds_no = 0; print_to_files && ds_no < max_blocks_to_print; ++ds_no)
{
auto ds = tp_result.GetPartition(ds_no);
std::string topVolumeBranchFileName = std::string("TopVolumeBranch_Rank_") +
std::to_string(static_cast<int>(rank)) + std::string("_Block_") +
std::to_string(static_cast<int>(ds_no)) + std::string(".txt");
std::ofstream topVolumeBranchStream(topVolumeBranchFileName.c_str());
auto topVolBranchGRId = ds.GetField("TopVolumeBranchGlobalRegularIds")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()
.ReadPortal();
auto topVolBranchVolume = ds.GetField("TopVolumeBranchVolume")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()
.ReadPortal();
auto topVolBranchSaddleEpsilon = ds.GetField("TopVolumeBranchSaddleEpsilon")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()
.ReadPortal();
auto topVolBranchSaddleIsoValue = ds.GetField("TopVolumeBranchSaddleIsoValue")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<ValueType>>()
.ReadPortal();
vtkm::Id nSelectedBranches = topVolBranchGRId.GetNumberOfValues();
for (vtkm::Id branch = 0; branch < nSelectedBranches; ++branch)
{
topVolumeBranchStream << std::setw(12) << topVolBranchGRId.Get(branch)
<< std::setw(14) << topVolBranchVolume.Get(branch)
<< std::setw(5) << topVolBranchSaddleEpsilon.Get(branch)
<< std::setw(14) << topVolBranchSaddleIsoValue.Get(branch)
<< std::endl;
}
std::string isoValuesFileName = std::string("IsoValues_Rank_") +
std::to_string(static_cast<int>(rank)) + std::string("_Block_") +
std::to_string(static_cast<int>(ds_no)) + std::string(".txt");
std::ofstream isoValuesStream(isoValuesFileName.c_str());
for (vtkm::Id branch = 0; branch < nSelectedBranches; ++branch)
{
isoValuesStream << (topVolBranchSaddleIsoValue.Get(branch) +
(eps * topVolBranchSaddleEpsilon.Get(branch)))
<< " ";
}
isoValuesStream << std::endl;
}
}
}
else
{

@ -12,17 +12,21 @@ set(scalar_topology_headers
ContourTreeUniformAugmented.h
ContourTreeUniformDistributed.h
DistributedBranchDecompositionFilter.h
SelectTopVolumeContoursFilter.h
)
set(scalar_topology_sources
internal/BranchDecompositionBlock.cxx
internal/SelectTopVolumeContoursBlock.cxx
internal/ComputeBlockIndices.cxx
internal/ComputeDistributedBranchDecompositionFunctor.cxx
internal/SelectTopVolumeContoursFunctor.cxx
internal/ExchangeBranchEndsFunctor.cxx
ContourTreeUniform.cxx
ContourTreeUniformAugmented.cxx
ContourTreeUniformDistributed.cxx
DistributedBranchDecompositionFilter.cxx
SelectTopVolumeContoursFilter.cxx
)
vtkm_library(

@ -832,6 +832,15 @@ VTKM_CONT void ContourTreeUniformDistributed::DoPostExecute(
"from information in CellSetStructured.");
diyBounds = vtkm::filter::scalar_topology::internal::ComputeBlockIndices(
input, diyDivisions, vtkmdiyLocalBlockGids);
// Set BlocksPerDimension fromn diyDivisions result as add them
// as information to the output data set for use in subsequent
// filters
this->BlocksPerDimension = vtkm::Id3{ 1, 1, 1 };
for (unsigned int d = 0; d < diyDivisions.size(); ++d)
{
this->BlocksPerDimension[d] = diyDivisions[d];
}
}
else
{

@ -331,7 +331,6 @@ VTKM_CONT vtkm::cont::PartitionedDataSet DistributedBranchDecompositionFilter::D
outputDataSets[ds_no] = input.GetPartition(ds_no);
}
branch_decomposition_master.foreach (
[&](BranchDecompositionBlock* b, const vtkmdiy::Master::ProxyWithLink&) {
vtkm::cont::Field branchRootField(
@ -347,6 +346,51 @@ VTKM_CONT vtkm::cont::PartitionedDataSet DistributedBranchDecompositionFilter::D
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.LowerEndGRId);
outputDataSets[b->LocalBlockNo].AddField(LowerEndGRIdField);
vtkm::cont::Field UpperEndIntrinsicVolume(
"UpperEndIntrinsicVolume",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.UpperEndIntrinsicVolume);
outputDataSets[b->LocalBlockNo].AddField(UpperEndIntrinsicVolume);
vtkm::cont::Field UpperEndDependentVolume(
"UpperEndDependentVolume",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.UpperEndDependentVolume);
outputDataSets[b->LocalBlockNo].AddField(UpperEndDependentVolume);
vtkm::cont::Field LowerEndIntrinsicVolume(
"LowerEndIntrinsicVolume",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.LowerEndIntrinsicVolume);
outputDataSets[b->LocalBlockNo].AddField(LowerEndIntrinsicVolume);
vtkm::cont::Field LowerEndDependentVolume(
"LowerEndDependentVolume",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.LowerEndDependentVolume);
outputDataSets[b->LocalBlockNo].AddField(LowerEndDependentVolume);
vtkm::cont::Field LowerEndSuperarcId("LowerEndSuperarcId",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.LowerEndSuperarcId);
outputDataSets[b->LocalBlockNo].AddField(LowerEndSuperarcId);
vtkm::cont::Field UpperEndSuperarcId("UpperEndSuperarcId",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.UpperEndSuperarcId);
outputDataSets[b->LocalBlockNo].AddField(UpperEndSuperarcId);
vtkm::cont::Field LowerEndValue("LowerEndValue",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.LowerEndValue);
outputDataSets[b->LocalBlockNo].AddField(LowerEndValue);
vtkm::cont::Field UpperEndValue("UpperEndValue",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.UpperEndValue);
outputDataSets[b->LocalBlockNo].AddField(UpperEndValue);
vtkm::cont::Field BranchRoot("BranchRoot",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.BranchRoot);
outputDataSets[b->LocalBlockNo].AddField(BranchRoot);
vtkm::cont::Field BranchRootGRId("BranchRootGRId",
vtkm::cont::Field::Association::WholeDataSet,
b->VolumetricBranchDecomposer.BranchRootGRId);
outputDataSets[b->LocalBlockNo].AddField(BranchRootGRId);
});
timingsStream << " " << std::setw(38) << std::left

@ -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 <vtkm/cont/EnvironmentTracker.h>
#include <vtkm/cont/ErrorFilterExecution.h>
#include <vtkm/filter/scalar_topology/SelectTopVolumeContoursFilter.h>
#include <vtkm/filter/scalar_topology/internal/SelectTopVolumeContoursBlock.h>
#include <vtkm/filter/scalar_topology/internal/SelectTopVolumeContoursFunctor.h>
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/ArrayTransforms.h>
// vtkm includes
#include <vtkm/cont/Timer.h>
// DIY includes
// clang-format off
VTKM_THIRDPARTY_PRE_INCLUDE
#include <vtkm/thirdparty/diy/Configure.h>
#include <vtkm/thirdparty/diy/diy.h>
VTKM_THIRDPARTY_POST_INCLUDE
// clang-format on
namespace vtkm
{
namespace filter
{
namespace scalar_topology
{
VTKM_CONT vtkm::cont::DataSet SelectTopVolumeContoursFilter::DoExecute(const vtkm::cont::DataSet&)
{
throw vtkm::cont::ErrorFilterExecution(
"SelectTopVolumeContoursFilter expects PartitionedDataSet as input.");
}
VTKM_CONT vtkm::cont::PartitionedDataSet SelectTopVolumeContoursFilter::DoExecutePartitions(
const vtkm::cont::PartitionedDataSet& input)
{
auto comm = vtkm::cont::EnvironmentTracker::GetCommunicator();
int rank = comm.rank();
int size = comm.size();
using SelectTopVolumeContoursBlock =
vtkm::filter::scalar_topology::internal::SelectTopVolumeContoursBlock;
vtkmdiy::Master branch_top_volume_master(comm,
1, // Use 1 thread, VTK-M will do the treading
-1, // All blocks in memory
0, // No create function
SelectTopVolumeContoursBlock::Destroy);
auto firstDS = input.GetPartition(0);
vtkm::Id3 firstPointDimensions, firstGlobalPointDimensions, firstGlobalPointIndexStart;
firstDS.GetCellSet().CastAndCallForTypes<VTKM_DEFAULT_CELL_SET_LIST_STRUCTURED>(
vtkm::worklet::contourtree_augmented::GetLocalAndGlobalPointDimensions(),
firstPointDimensions,
firstGlobalPointDimensions,
firstGlobalPointIndexStart);
int numDims = firstGlobalPointDimensions[2] > 1 ? 3 : 2;
auto vtkmBlocksPerDimensionRP = input.GetPartition(0)
.GetField("vtkmBlocksPerDimension")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()
.ReadPortal();
int globalNumberOfBlocks = 1;
for (vtkm::IdComponent d = 0; d < static_cast<vtkm::IdComponent>(numDims); ++d)
{
globalNumberOfBlocks *= static_cast<int>(vtkmBlocksPerDimensionRP.Get(d));
}
vtkmdiy::DynamicAssigner assigner(comm, size, globalNumberOfBlocks);
for (vtkm::Id localBlockIndex = 0; localBlockIndex < input.GetNumberOfPartitions();
++localBlockIndex)
{
const vtkm::cont::DataSet& ds = input.GetPartition(localBlockIndex);
int globalBlockId = static_cast<int>(
vtkm::cont::ArrayGetValue(0,
ds.GetField("vtkmGlobalBlockId")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()));
SelectTopVolumeContoursBlock* b =
new SelectTopVolumeContoursBlock(localBlockIndex, globalBlockId);
branch_top_volume_master.add(globalBlockId, b, new vtkmdiy::Link());
assigner.set_rank(rank, globalBlockId);
}
vtkmdiy::fix_links(branch_top_volume_master, assigner);
branch_top_volume_master.foreach (
[&](SelectTopVolumeContoursBlock* b, const vtkmdiy::Master::ProxyWithLink&) {
const auto& globalSize = firstGlobalPointDimensions;
vtkm::Id totalVolume = globalSize[0] * globalSize[1] * globalSize[2];
const vtkm::cont::DataSet& ds = input.GetPartition(b->LocalBlockNo);
b->SortBranchByVolume(ds, totalVolume);
// copy the top volume branches into a smaller array
// we skip index 0 because it must be the main branch (which has the highest volume)
vtkm::Id nActualSavedBranches =
std::min(this->nSavedBranches, b->SortedBranchByVolume.GetNumberOfValues() - 1);
vtkm::worklet::contourtree_augmented::IdArrayType topVolumeBranch;
vtkm::cont::Algorithm::CopySubRange(
b->SortedBranchByVolume, 1, nActualSavedBranches, topVolumeBranch);
auto branchRootGRId =
ds.GetField("BranchRootGRId").GetData().AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
branchRootGRId, topVolumeBranch, b->TopVolumeBranchRootGRId);
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
branchRootGRId, topVolumeBranch, b->TopVolumeBranchRootGRId);
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
b->BranchVolume, topVolumeBranch, b->TopVolumeBranchVolume);
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
b->BranchSaddleEpsilon, topVolumeBranch, b->TopVolumeBranchSaddleEpsilon);
auto resolveArray = [&](const auto& inArray) {
using InArrayHandleType = std::decay_t<decltype(inArray)>;
InArrayHandleType topVolBranchSaddleIsoValue;
vtkm::worklet::contourtree_augmented::PermuteArrayWithRawIndex<InArrayHandleType>(
inArray, topVolumeBranch, topVolBranchSaddleIsoValue);
b->TopVolumeBranchSaddleIsoValue = topVolBranchSaddleIsoValue;
};
b->BranchSaddleIsoValue
.CastAndCallForTypes<vtkm::TypeListScalarAll, vtkm::cont::StorageListBasic>(resolveArray);
});
// We apply all-to-all broadcast to collect the top nSavedBranches branches by volume
vtkmdiy::all_to_all(
branch_top_volume_master,
assigner,
vtkm::filter::scalar_topology::internal::SelectTopVolumeContoursFunctor(this->nSavedBranches));
// For each block, we compute the get the extracted isosurface for every selected branch
// storing format: key (branch ID) - Value (list of meshes in the isosurface)
std::vector<vtkm::cont::DataSet> outputDataSets(input.GetNumberOfPartitions());
branch_top_volume_master.foreach (
[&](SelectTopVolumeContoursBlock* b, const vtkmdiy::Master::ProxyWithLink&) {
vtkm::cont::Field TopVolBranchGRIdField("TopVolumeBranchGlobalRegularIds",
vtkm::cont::Field::Association::WholeDataSet,
b->TopVolumeBranchRootGRId);
outputDataSets[b->LocalBlockNo].AddField(TopVolBranchGRIdField);
vtkm::cont::Field TopVolBranchVolumeField("TopVolumeBranchVolume",
vtkm::cont::Field::Association::WholeDataSet,
b->TopVolumeBranchVolume);
outputDataSets[b->LocalBlockNo].AddField(TopVolBranchVolumeField);
vtkm::cont::Field TopVolBranchSaddleEpsilonField("TopVolumeBranchSaddleEpsilon",
vtkm::cont::Field::Association::WholeDataSet,
b->TopVolumeBranchSaddleEpsilon);
outputDataSets[b->LocalBlockNo].AddField(TopVolBranchSaddleEpsilonField);
auto resolveArray = [&](const auto& inArray) {
vtkm::cont::Field TopVolBranchSaddleIsoValueField(
"TopVolumeBranchSaddleIsoValue", vtkm::cont::Field::Association::WholeDataSet, inArray);
outputDataSets[b->LocalBlockNo].AddField(TopVolBranchSaddleIsoValueField);
};
b->TopVolumeBranchSaddleIsoValue
.CastAndCallForTypes<vtkm::TypeListScalarAll, vtkm::cont::StorageListBasic>(resolveArray);
});
return vtkm::cont::PartitionedDataSet{ outputDataSets };
}
} // namespace scalar_topology
} // namespace filter
} // namespace vtkm

@ -0,0 +1,83 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
#ifndef vtk_m_filter_scalar_topology_SelectTopVolumeContoursFilter_h
#define vtk_m_filter_scalar_topology_SelectTopVolumeContoursFilter_h
#include <vtkm/filter/Filter.h>
#include <vtkm/filter/scalar_topology/vtkm_filter_scalar_topology_export.h>
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/Types.h>
namespace vtkm
{
namespace filter
{
namespace scalar_topology
{
/// \brief Compute branch decompostion from distributed contour tree
class VTKM_FILTER_SCALAR_TOPOLOGY_EXPORT SelectTopVolumeContoursFilter : public vtkm::filter::Filter
{
public:
VTKM_CONT SelectTopVolumeContoursFilter() = default;
VTKM_CONT void SetSavedBranches(const vtkm::Id& numBranches)
{
this->nSavedBranches = numBranches;
}
private:
VTKM_CONT vtkm::cont::DataSet DoExecute(const vtkm::cont::DataSet&) override;
VTKM_CONT vtkm::cont::PartitionedDataSet DoExecutePartitions(
const vtkm::cont::PartitionedDataSet& inData) override;
vtkm::Id nSavedBranches;
vtkm::worklet::contourtree_augmented::IdArrayType BranchVolume;
vtkm::worklet::contourtree_augmented::IdArrayType BranchSaddleEpsilon;
vtkm::worklet::contourtree_augmented::IdArrayType SortedBranchByVolume;
vtkm::cont::UnknownArrayHandle BranchSaddleIsoValue;
};
} // namespace scalar_topology
} // namespace worklet
} // namespace vtkm
#endif

@ -10,8 +10,10 @@
set(headers
BranchDecompositionBlock.h
SelectTopVolumeContoursBlock.h
ComputeBlockIndices.h
ComputeDistributedBranchDecompositionFunctor.h
SelectTopVolumeContoursFunctor.h
ExchangeBranchEndsFunctor.h
)
#-----------------------------------------------------------------------------

@ -0,0 +1,229 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#include <vtkm/cont/ArrayCopy.h>
#include <vtkm/filter/scalar_topology/internal/SelectTopVolumeContoursBlock.h>
#include <vtkm/filter/scalar_topology/worklet/select_top_volume_contours/ClarifyBranchEndSupernodeTypeWorklet.h>
#include <vtkm/filter/scalar_topology/worklet/select_top_volume_contours/GetBranchVolumeWorklet.h>
#include <vtkm/filter/scalar_topology/worklet/select_top_volume_contours/UpdateInfoByBranchDirectionWorklet.h>
#ifdef DEBUG_PRINT
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/PrintVectors.h>
#endif
namespace vtkm
{
namespace filter
{
namespace scalar_topology
{
namespace internal
{
SelectTopVolumeContoursBlock::SelectTopVolumeContoursBlock(vtkm::Id localBlockNo, int globalBlockId)
: LocalBlockNo(localBlockNo)
, GlobalBlockId(globalBlockId)
{
}
void SelectTopVolumeContoursBlock::SortBranchByVolume(
const vtkm::cont::DataSet& hierarchicalTreeDataSet,
const vtkm::Id totalVolume)
{
/// Pipeline to compute the branch volume
/// 1. check both ends of the branch. If both leaves, then main branch, volume = totalVolume
/// 2. for other branches, check the direction of the inner superarc
/// branch volume = (inner superarc points to the senior-most node) ?
/// dependentVolume[innerSuperarc] :
/// reverseVolume[innerSuperarc]
/// NOTE: reverseVolume = totalVolume - dependentVolume + intrinsicVolume
// Generally, if ending superarc has intrinsicVol == dependentVol, then it is a leaf node
vtkm::cont::ArrayHandle<bool> isLowerLeaf;
vtkm::cont::ArrayHandle<bool> isUpperLeaf;
auto upperEndIntrinsicVolume = hierarchicalTreeDataSet.GetField("UpperEndIntrinsicVolume")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
auto upperEndDependentVolume = hierarchicalTreeDataSet.GetField("UpperEndDependentVolume")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
auto lowerEndIntrinsicVolume = hierarchicalTreeDataSet.GetField("LowerEndIntrinsicVolume")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
auto lowerEndDependentVolume = hierarchicalTreeDataSet.GetField("LowerEndDependentVolume")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
auto lowerEndSuperarcId = hierarchicalTreeDataSet.GetField("LowerEndSuperarcId")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
auto upperEndSuperarcId = hierarchicalTreeDataSet.GetField("UpperEndSuperarcId")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
auto branchRoot = hierarchicalTreeDataSet.GetField("BranchRoot")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>();
vtkm::cont::Algorithm::Transform(
upperEndIntrinsicVolume, upperEndDependentVolume, isUpperLeaf, vtkm::Equal());
vtkm::cont::Algorithm::Transform(
lowerEndIntrinsicVolume, lowerEndDependentVolume, isLowerLeaf, vtkm::Equal());
// NOTE: special cases (one-superarc branches) exist
// if the upper end superarc == lower end superarc == branch root superarc
// then it's probably not a leaf-leaf branch (Both equality has to be satisfied!)
// exception: the entire domain has only one superarc (intrinsic == dependent == total - 1)
// then it is a leaf-leaf branch
vtkm::cont::Invoker invoke;
vtkm::worklet::scalar_topology::select_top_volume_contours::ClarifyBranchEndSupernodeTypeWorklet
clarifyNodeTypeWorklet(totalVolume);
invoke(clarifyNodeTypeWorklet,
lowerEndSuperarcId,
lowerEndIntrinsicVolume,
upperEndSuperarcId,
upperEndIntrinsicVolume,
branchRoot,
isLowerLeaf,
isUpperLeaf);
vtkm::cont::UnknownArrayHandle upperEndValue =
hierarchicalTreeDataSet.GetField("UpperEndValue").GetData();
// Based on the direction info of the branch, store epsilon direction and isovalue of the saddle
auto resolveArray = [&](const auto& inArray) {
using InArrayHandleType = std::decay_t<decltype(inArray)>;
using ValueType = typename InArrayHandleType::ValueType;
vtkm::cont::ArrayHandle<ValueType> branchSaddleIsoValue;
branchSaddleIsoValue.Allocate(isLowerLeaf.GetNumberOfValues());
this->BranchSaddleEpsilon.Allocate(isLowerLeaf.GetNumberOfValues());
vtkm::worklet::scalar_topology::select_top_volume_contours::UpdateInfoByBranchDirectionWorklet<
ValueType>
updateInfoWorklet;
auto lowerEndValue = hierarchicalTreeDataSet.GetField("LowerEndValue")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<ValueType>>();
invoke(updateInfoWorklet,
isLowerLeaf,
isUpperLeaf,
inArray,
lowerEndValue,
this->BranchSaddleEpsilon,
branchSaddleIsoValue);
this->BranchSaddleIsoValue = branchSaddleIsoValue;
};
upperEndValue.CastAndCallForTypes<vtkm::TypeListScalarAll, vtkm::cont::StorageListBasic>(
resolveArray);
vtkm::worklet::contourtree_augmented::IdArrayType branchVolume;
vtkm::worklet::scalar_topology::select_top_volume_contours::GetBranchVolumeWorklet
getBranchVolumeWorklet(totalVolume);
invoke(getBranchVolumeWorklet, // worklet
lowerEndSuperarcId, // input
lowerEndIntrinsicVolume, // input
lowerEndDependentVolume, // input
upperEndSuperarcId, // input
upperEndIntrinsicVolume, // input
upperEndDependentVolume, // input
isLowerLeaf,
isUpperLeaf,
branchVolume); // output
#ifdef DEBUG_PRINT
std::stringstream resultStream;
resultStream << "Branch Volume In The Block" << std::endl;
const vtkm::Id nVolume = branchVolume.GetNumberOfValues();
vtkm::worklet::contourtree_augmented::PrintHeader(nVolume, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"BranchVolume", branchVolume, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices("isLowerLeaf", isLowerLeaf, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices("isUpperLeaf", isUpperLeaf, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"LowerEndIntrinsicVol", lowerEndIntrinsicVolume, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"LowerEndDependentVol", lowerEndDependentVolume, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"UpperEndIntrinsicVol", upperEndIntrinsicVolume, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"UpperEndDependentVol", upperEndDependentVolume, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"LowerEndSuperarc", lowerEndSuperarcId, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"UpperEndSuperarc", upperEndSuperarcId, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices("BranchRoot", branchRoot, -1, resultStream);
resultStream << std::endl;
VTKM_LOG_S(vtkm::cont::LogLevel::Info, resultStream.str());
#endif
vtkm::cont::Algorithm::Copy(branchVolume, this->BranchVolume);
const vtkm::Id nBranches = lowerEndSuperarcId.GetNumberOfValues();
vtkm::cont::ArrayHandleIndex branchesIdx(nBranches);
vtkm::worklet::contourtree_augmented::IdArrayType sortedBranches;
vtkm::cont::Algorithm::Copy(branchesIdx, sortedBranches);
// sort the branch volume
vtkm::cont::Algorithm::SortByKey(branchVolume, sortedBranches, vtkm::SortGreater());
vtkm::cont::Algorithm::Copy(sortedBranches, this->SortedBranchByVolume);
}
} // namespace internal
} // namespace scalar_topology
} // namespace filter
} // namespace vtkm

@ -0,0 +1,98 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#ifndef vtk_m_filter_scalar_topology_internal_SelectTopVolumeContoursBlock_h
#define vtk_m_filter_scalar_topology_internal_SelectTopVolumeContoursBlock_h
#include <vtkm/cont/DataSet.h>
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/Types.h>
namespace vtkm
{
namespace filter
{
namespace scalar_topology
{
namespace internal
{
struct SelectTopVolumeContoursBlock
{
SelectTopVolumeContoursBlock(vtkm::Id localBlockNo, int globalBlockId);
void SortBranchByVolume(const vtkm::cont::DataSet& hierarchicalTreeDataSet,
const vtkm::Id totalVolume);
// Block metadata
vtkm::Id LocalBlockNo;
int GlobalBlockId; // TODO/FIXME: Check whether really needed. Possibly only during debugging
vtkm::worklet::contourtree_augmented::IdArrayType BranchVolume;
vtkm::worklet::contourtree_augmented::IdArrayType BranchSaddleEpsilon;
vtkm::worklet::contourtree_augmented::IdArrayType SortedBranchByVolume;
vtkm::cont::UnknownArrayHandle BranchSaddleIsoValue;
// Output Datasets.
vtkm::worklet::contourtree_augmented::IdArrayType TopVolumeBranchRootGRId;
vtkm::worklet::contourtree_augmented::IdArrayType TopVolumeBranchVolume;
vtkm::cont::UnknownArrayHandle TopVolumeBranchSaddleIsoValue;
vtkm::worklet::contourtree_augmented::IdArrayType TopVolumeBranchSaddleEpsilon;
// Destroy function allowing DIY to own blocks and clean them up after use
static void Destroy(void* b) { delete static_cast<SelectTopVolumeContoursBlock*>(b); }
};
} // namespace internal
} // namespace scalar_topology
} // namespace filter
} // namespace vtkm
#endif

@ -0,0 +1,384 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#include <vtkm/filter/scalar_topology/internal/SelectTopVolumeContoursFunctor.h>
#include <vtkm/filter/scalar_topology/worklet/branch_decomposition/hierarchical_volumetric_branch_decomposer/GetOuterEndWorklet.h>
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/ArrayTransforms.h>
#include <vtkm/filter/scalar_topology/worklet/select_top_volume_contours/BranchVolumeComparator.h>
#include <vtkm/Types.h>
#ifdef DEBUG_PRINT
#define DEBUG_PRINT_COMBINED_HIGH_VOLUME_BRANCH
#include <vtkm/filter/scalar_topology/worklet/contourtree_augmented/PrintVectors.h>
#endif
namespace vtkm
{
namespace filter
{
namespace scalar_topology
{
namespace internal
{
void SelectTopVolumeContoursFunctor::operator()(
SelectTopVolumeContoursBlock* b,
const vtkmdiy::ReduceProxy& rp // communication proxy
) const
{
if (this->nSavedBranches < 1)
return;
// Get our rank and DIY id
//const vtkm::Id rank = vtkm::cont::EnvironmentTracker::GetCommunicator().rank();
const auto selfid = rp.gid();
// Aliases to reduce verbosity
using IdArrayType = vtkm::worklet::contourtree_augmented::IdArrayType;
vtkm::cont::Invoker invoke;
if (rp.in_link().size() == 0)
{
for (int cc = 0; cc < rp.out_link().size(); ++cc)
{
auto target = rp.out_link().target(cc);
if (target.gid != selfid)
{
#ifdef DEBUG_PRINT_COMBINED_HIGH_VOLUME_BRANCH
rp.enqueue(target, b->GlobalBlockId);
VTKM_LOG_S(vtkm::cont::LogLevel::Info,
"Block " << b->GlobalBlockId << " enqueue to Block " << target.gid);
#endif
auto topVolBranchRootGRIdPortal = b->TopVolumeBranchRootGRId.ReadPortal();
auto topVolBranchVolumePortal = b->TopVolumeBranchVolume.ReadPortal();
auto topVolBranchSaddleEpsilonPortal = b->TopVolumeBranchSaddleEpsilon.ReadPortal();
vtkm::Id nBranches = topVolBranchRootGRIdPortal.GetNumberOfValues();
rp.enqueue(target, nBranches);
for (vtkm::Id branch = 0; branch < nBranches; ++branch)
rp.enqueue(target, topVolBranchRootGRIdPortal.Get(branch));
for (vtkm::Id branch = 0; branch < nBranches; ++branch)
rp.enqueue(target, topVolBranchVolumePortal.Get(branch));
for (vtkm::Id branch = 0; branch < nBranches; ++branch)
rp.enqueue(target, topVolBranchSaddleEpsilonPortal.Get(branch));
auto resolveArray = [&](const auto& inArray) {
using InArrayHandleType = std::decay_t<decltype(inArray)>;
using ValueType = typename InArrayHandleType::ValueType;
auto topVolBranchSaddleIsoValuePortal = inArray.ReadPortal();
for (vtkm::Id branch = 0; branch < nBranches; ++branch)
rp.enqueue<ValueType>(target, topVolBranchSaddleIsoValuePortal.Get(branch));
};
b->TopVolumeBranchSaddleIsoValue
.CastAndCallForTypes<vtkm::TypeListScalarAll, vtkm::cont::StorageListBasic>(resolveArray);
// rp.enqueue(target, b->TopVolumeBranchRootGRId);
// rp.enqueue(target, b->TopVolumeBranchVolume);
}
}
}
else
{
for (int i = 0; i < rp.in_link().size(); ++i)
{
int ingid = rp.in_link().target(i).gid;
if (ingid == selfid)
continue;
// copy incoming to the block
#ifdef DEBUG_PRINT_COMBINED_HIGH_VOLUME_BRANCH
int incomingGlobalBlockId;
rp.dequeue(ingid, incomingGlobalBlockId);
VTKM_LOG_S(vtkm::cont::LogLevel::Info,
"Combining local block " << b->GlobalBlockId << " with incoming block "
<< incomingGlobalBlockId);
#endif
// dequeue the data from other blocks.
// nIncomingBranches
// array of incoming branch global regular ID
// array of incoming branch volume
// array of branch epsilon direction
// array of branch saddle end value
// TODO/FIXME: This is a workaround for a bug in DIY/vtk-m.
// Replace with dequeuing ArrayHandles once bug is fixed.
vtkm::Id nIncoming;
rp.dequeue(ingid, nIncoming);
IdArrayType incomingTopVolBranchGRId;
incomingTopVolBranchGRId.Allocate(nIncoming);
auto incomingTopVolBranchGRIdPortal = incomingTopVolBranchGRId.WritePortal();
for (vtkm::Id branch = 0; branch < nIncoming; ++branch)
{
vtkm::Id incomingTmpBranchGRId;
rp.dequeue(ingid, incomingTmpBranchGRId);
incomingTopVolBranchGRIdPortal.Set(branch, incomingTmpBranchGRId);
}
// TODO/FIXME: This is a workaround for a bug in DIY/vtk-m.
// Replace with dequeuing ArrayHandles once bug is fixed.
// rp.dequeue(ingid, incomingTopVolBranchGRId);
IdArrayType incomingTopVolBranchVolume;
incomingTopVolBranchVolume.Allocate(nIncoming);
auto incomingTopVolBranchVolumePortal = incomingTopVolBranchVolume.WritePortal();
for (vtkm::Id branch = 0; branch < nIncoming; ++branch)
{
vtkm::Id incomingTmpBranchVolume;
rp.dequeue(ingid, incomingTmpBranchVolume);
incomingTopVolBranchVolumePortal.Set(branch, incomingTmpBranchVolume);
}
// TODO/FIXME: This is a workaround for a bug in DIY/vtk-m.
// Replace with dequeuing ArrayHandles once bug is fixed.
// rp.dequeue(ingid, incomingTopVolBranchVolume);
IdArrayType incomingTopVolBranchSaddleEpsilon;
incomingTopVolBranchSaddleEpsilon.Allocate(nIncoming);
auto incomingTopVolBranchSaddleEpsilonPortal =
incomingTopVolBranchSaddleEpsilon.WritePortal();
for (vtkm::Id branch = 0; branch < nIncoming; ++branch)
{
vtkm::Id incomingTmpBranchSaddleEpsilon;
rp.dequeue(ingid, incomingTmpBranchSaddleEpsilon);
incomingTopVolBranchSaddleEpsilonPortal.Set(branch, incomingTmpBranchSaddleEpsilon);
}
// TODO/FIXME: This is a workaround for a bug in DIY/vtk-m.
// Replace with dequeuing ArrayHandles once bug is fixed.
// rp.dequeue(ingid, incomingTopVolBranchSaddleEpsilon);
auto resolveArray = [&](auto& inArray) {
using InArrayHandleType = std::decay_t<decltype(inArray)>;
using ValueType = typename InArrayHandleType::ValueType;
InArrayHandleType incomingTopVolBranchSaddleIsoValue;
incomingTopVolBranchSaddleIsoValue.Allocate(nIncoming);
auto incomingTopVolBranchSaddleIsoValuePortal =
incomingTopVolBranchSaddleIsoValue.WritePortal();
for (vtkm::Id branch = 0; branch < nIncoming; ++branch)
{
ValueType incomingSaddleValue;
rp.dequeue<ValueType>(ingid, incomingSaddleValue);
incomingTopVolBranchSaddleIsoValuePortal.Set(branch, incomingSaddleValue);
}
// TODO/FIXME: This is a workaround for a bug in DIY/vtk-m.
// Replace with dequeuing ArrayHandles once bug is fixed.
// rp.dequeue<InArrayHandleType>(ingid, incomingTopVolBranchSaddleIsoValue);
vtkm::Id nSelf = b->TopVolumeBranchRootGRId.GetNumberOfValues();
#ifdef DEBUG_PRINT_COMBINED_HIGH_VOLUME_BRANCH
VTKM_LOG_S(vtkm::cont::LogLevel::Info,
"nIncoming = " << nIncoming << ", nSelf = " << nSelf);
{
std::stringstream rs;
vtkm::worklet::contourtree_augmented::PrintHeader(nIncoming, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"incomingTopBranchId", incomingTopVolBranchGRId, -1, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"incomingTopBranchVol", incomingTopVolBranchVolume, -1, rs);
vtkm::worklet::contourtree_augmented::PrintValues<ValueType>(
"incomingSaddleVal", incomingTopVolBranchSaddleIsoValue, -1, rs);
vtkm::worklet::contourtree_augmented::PrintHeader(nSelf, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"selfTopBranchId", b->TopVolumeBranchRootGRId, -1, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"selfTopBranchVol", b->TopVolumeBranchVolume, -1, rs);
vtkm::worklet::contourtree_augmented::PrintValues<ValueType>(
"selfTopSaddleVal", inArray, -1, rs);
VTKM_LOG_S(vtkm::cont::LogLevel::Info, rs.str());
}
#endif
// merge incoming branches with self branches
IdArrayType mergedTopVolBranchGRId;
IdArrayType mergedTopVolBranchVolume;
IdArrayType mergedTopVolBranchSaddleEpsilon;
InArrayHandleType mergedTopVolBranchSaddleIsoValue;
mergedTopVolBranchGRId.Allocate(nIncoming + nSelf);
mergedTopVolBranchVolume.Allocate(nIncoming + nSelf);
mergedTopVolBranchSaddleEpsilon.Allocate(nIncoming + nSelf);
mergedTopVolBranchSaddleIsoValue.Allocate(nIncoming + nSelf);
vtkm::cont::Algorithm::CopySubRange(
incomingTopVolBranchGRId, 0, nIncoming, mergedTopVolBranchGRId, 0);
vtkm::cont::Algorithm::CopySubRange(
incomingTopVolBranchVolume, 0, nIncoming, mergedTopVolBranchVolume, 0);
vtkm::cont::Algorithm::CopySubRange(
incomingTopVolBranchSaddleEpsilon, 0, nIncoming, mergedTopVolBranchSaddleEpsilon, 0);
vtkm::cont::Algorithm::CopySubRange<ValueType, ValueType>(
incomingTopVolBranchSaddleIsoValue, 0, nIncoming, mergedTopVolBranchSaddleIsoValue, 0);
vtkm::cont::Algorithm::CopySubRange(
b->TopVolumeBranchRootGRId, 0, nSelf, mergedTopVolBranchGRId, nIncoming);
vtkm::cont::Algorithm::CopySubRange(
b->TopVolumeBranchVolume, 0, nSelf, mergedTopVolBranchVolume, nIncoming);
vtkm::cont::Algorithm::CopySubRange(
b->TopVolumeBranchSaddleEpsilon, 0, nSelf, mergedTopVolBranchSaddleEpsilon, nIncoming);
vtkm::cont::Algorithm::CopySubRange<ValueType, ValueType>(
inArray, 0, nSelf, mergedTopVolBranchSaddleIsoValue, nIncoming);
// Sort all branches (incoming + self) based on volume
// sorting key: (volume, branch global regular ID)
// the highest volume comes first, the lowest branch GR ID comes first
vtkm::cont::ArrayHandleIndex mergedBranchId(nIncoming + nSelf);
IdArrayType sortedBranchId;
vtkm::cont::Algorithm::Copy(mergedBranchId, sortedBranchId);
vtkm::worklet::scalar_topology::select_top_volume_contours::BranchVolumeComparator
branchVolumeComparator(mergedTopVolBranchGRId, mergedTopVolBranchVolume);
vtkm::cont::Algorithm::Sort(sortedBranchId, branchVolumeComparator);
// permute the branch information based on sorting
IdArrayType permutedTopVolBranchGRId;
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id, IdArrayType>(
mergedTopVolBranchGRId, sortedBranchId, permutedTopVolBranchGRId);
IdArrayType permutedTopVolBranchVolume;
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id, IdArrayType>(
mergedTopVolBranchVolume, sortedBranchId, permutedTopVolBranchVolume);
IdArrayType permutedTopVolBranchSaddleEpsilon;
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id, IdArrayType>(
mergedTopVolBranchSaddleEpsilon, sortedBranchId, permutedTopVolBranchSaddleEpsilon);
InArrayHandleType permutedTopVolBranchSaddleIsoValue;
vtkm::worklet::contourtree_augmented::PermuteArrayWithRawIndex<InArrayHandleType>(
mergedTopVolBranchSaddleIsoValue, sortedBranchId, permutedTopVolBranchSaddleIsoValue);
#ifdef DEBUG_PRINT_COMBINED_HIGH_VOLUME_BRANCH
{
std::stringstream rs;
vtkm::worklet::contourtree_augmented::PrintHeader(nIncoming + nSelf, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"permutedTopBranchId", permutedTopVolBranchGRId, -1, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"permutedTopBranchVol", permutedTopVolBranchVolume, -1, rs);
vtkm::worklet::contourtree_augmented::PrintValues<ValueType>(
"permutedTopSaddleVal", permutedTopVolBranchSaddleIsoValue, -1, rs);
}
#endif
// there may be duplicate branches. We remove duplicate branches based on global regular IDs
// We can reuse the filter from removing duplicate branches in the process of collecting branches
IdArrayType oneIfUniqueBranch;
oneIfUniqueBranch.Allocate(nIncoming + nSelf);
vtkm::worklet::scalar_topology::hierarchical_volumetric_branch_decomposer::
OneIfBranchEndWorklet oneIfUniqueWorklet;
invoke(oneIfUniqueWorklet, mergedBranchId, permutedTopVolBranchGRId, oneIfUniqueBranch);
// Remove duplicate
IdArrayType mergedUniqueBranchGRId;
IdArrayType mergedUniqueBranchVolume;
IdArrayType mergedUniqueBranchSaddleEpsilon;
InArrayHandleType mergedUniqueBranchSaddleIsoValue;
vtkm::cont::Algorithm::CopyIf(
permutedTopVolBranchGRId, oneIfUniqueBranch, mergedUniqueBranchGRId);
vtkm::cont::Algorithm::CopyIf(
permutedTopVolBranchVolume, oneIfUniqueBranch, mergedUniqueBranchVolume);
vtkm::cont::Algorithm::CopyIf(
permutedTopVolBranchSaddleEpsilon, oneIfUniqueBranch, mergedUniqueBranchSaddleEpsilon);
vtkm::cont::Algorithm::CopyIf(
permutedTopVolBranchSaddleIsoValue, oneIfUniqueBranch, mergedUniqueBranchSaddleIsoValue);
vtkm::Id nMergedUnique = mergedUniqueBranchGRId.GetNumberOfValues();
#ifdef DEBUG_PRINT_COMBINED_HIGH_VOLUME_BRANCH
{
std::stringstream rs;
vtkm::worklet::contourtree_augmented::PrintHeader(nMergedUnique, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"mergedUniqueBranchId", mergedUniqueBranchGRId, -1, rs);
vtkm::worklet::contourtree_augmented::PrintIndices(
"mergedUniqueBranchVol", mergedUniqueBranchVolume, -1, rs);
vtkm::worklet::contourtree_augmented::PrintValues<ValueType>(
"mergedUniqueSaddleVal", mergedUniqueBranchSaddleIsoValue, -1, rs);
}
#endif
// After removing duplicate, if there are more branches than we need
// We only save the top nSavedBranches branches
if (nMergedUnique > this->nSavedBranches)
{
vtkm::cont::Algorithm::CopySubRange(
mergedUniqueBranchGRId, 0, this->nSavedBranches, b->TopVolumeBranchRootGRId);
vtkm::cont::Algorithm::CopySubRange(
mergedUniqueBranchVolume, 0, this->nSavedBranches, b->TopVolumeBranchVolume);
vtkm::cont::Algorithm::CopySubRange(mergedUniqueBranchSaddleEpsilon,
0,
this->nSavedBranches,
b->TopVolumeBranchSaddleEpsilon);
// InArrayHandleType subRangeUniqueBranchSaddleIsoValue;
inArray.Allocate(this->nSavedBranches);
vtkm::cont::Algorithm::CopySubRange(
mergedUniqueBranchSaddleIsoValue, 0, this->nSavedBranches, inArray);
// inArray = subRangeUniqueBranchSaddleIsoValue;
}
else
{
vtkm::cont::Algorithm::Copy(mergedUniqueBranchGRId, b->TopVolumeBranchRootGRId);
vtkm::cont::Algorithm::Copy(mergedUniqueBranchVolume, b->TopVolumeBranchVolume);
vtkm::cont::Algorithm::Copy(mergedUniqueBranchSaddleEpsilon,
b->TopVolumeBranchSaddleEpsilon);
inArray.Allocate(nMergedUnique);
vtkm::cont::Algorithm::Copy(mergedUniqueBranchSaddleIsoValue, inArray);
}
};
b->TopVolumeBranchSaddleIsoValue
.CastAndCallForTypes<vtkm::TypeListScalarAll, vtkm::cont::StorageListBasic>(resolveArray);
}
}
}
} // namespace internal
} // namespace scalar_topology
} // namespace filter
} // namespace vtkm

@ -0,0 +1,93 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#ifndef vtk_m_filter_scalar_topology_internal_SelectTopVolumeContoursFunctor_h
#define vtk_m_filter_scalar_topology_internal_SelectTopVolumeContoursFunctor_h
#include <vtkm/filter/scalar_topology/internal/SelectTopVolumeContoursBlock.h>
// clang-format off
VTKM_THIRDPARTY_PRE_INCLUDE
#include <vtkm/thirdparty/diy/diy.h>
VTKM_THIRDPARTY_POST_INCLUDE
// clang-format on
namespace vtkm
{
namespace filter
{
namespace scalar_topology
{
namespace internal
{
struct SelectTopVolumeContoursFunctor
{
SelectTopVolumeContoursFunctor(const vtkm::Id& nSB)
: nSavedBranches(nSB)
{
}
void operator()(SelectTopVolumeContoursBlock* b,
const vtkmdiy::ReduceProxy& rp // communication proxy
) const;
const vtkm::Id nSavedBranches;
};
} // namespace internal
} // namespace scalar_topology
} // namespace filter
} // namespace vtkm
#endif

@ -63,6 +63,7 @@
#include <vtkm/filter/MapFieldPermutation.h>
#include <vtkm/filter/scalar_topology/ContourTreeUniformDistributed.h>
#include <vtkm/filter/scalar_topology/DistributedBranchDecompositionFilter.h>
#include <vtkm/filter/scalar_topology/SelectTopVolumeContoursFilter.h>
#include <vtkm/filter/scalar_topology/testing/SuperArcHelper.h>
#include <vtkm/filter/scalar_topology/testing/VolumeHelper.h>
#include <vtkm/filter/scalar_topology/worklet/branch_decomposition/HierarchicalVolumetricBranchDecomposer.h>
@ -466,6 +467,14 @@ inline void TestContourTreeUniformDistributedBranchDecomposition8x9(int nBlocks,
augmentHierarchicalTree,
computeHierarchicalVolumetricBranchDecomposition);
using vtkm::filter::scalar_topology::SelectTopVolumeContoursFilter;
vtkm::Id numBranches = 2;
SelectTopVolumeContoursFilter tp_filter;
tp_filter.SetSavedBranches(numBranches);
auto tp_result = tp_filter.Execute(result);
if (vtkm::cont::EnvironmentTracker::GetCommunicator().rank() == 0)
{
@ -503,27 +512,99 @@ inline void TestContourTreeUniformDistributedBranchDecomposition8x9(int nBlocks,
std::sort(computed.begin(), computed.end());
std::sort(expected.begin(), expected.end());
if (computed == expected)
if (computed != expected)
{
std::cout << "Branch Decomposition: Results Match!" << std::endl;
return;
std::cout << "Branch Decomposition Results:" << std::endl;
std::cout << "Computed Contour Tree" << std::endl;
for (std::size_t i = 0; i < computed.size(); i++)
{
std::cout << std::setw(12) << computed[i].low << std::setw(14) << computed[i].high
<< std::endl;
}
std::cout << "Expected Contour Tree" << std::endl;
for (std::size_t i = 0; i < expected.size(); i++)
{
std::cout << std::setw(12) << expected[i].low << std::setw(14) << expected[i].high
<< std::endl;
}
VTKM_TEST_FAIL("Branch Decomposition Failed!");
}
std::cout << "Branch Decomposition Results:" << std::endl;
std::cout << "Computed Contour Tree" << std::endl;
for (std::size_t i = 0; i < computed.size(); i++)
std::cout << "Branch Decomposition: Results Match!" << std::endl;
for (vtkm::Id ds_no = 0; ds_no < result.GetNumberOfPartitions(); ++ds_no)
{
std::cout << std::setw(12) << computed[i].low << std::setw(14) << computed[i].high
<< std::endl;
auto ds = tp_result.GetPartition(ds_no);
auto topVolBranchGRId = ds.GetField("TopVolumeBranchGlobalRegularIds")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()
.ReadPortal();
auto topVolBranchVolume = ds.GetField("TopVolumeBranchVolume")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()
.ReadPortal();
auto topVolBranchSaddleEpsilon = ds.GetField("TopVolumeBranchSaddleEpsilon")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Id>>()
.ReadPortal();
auto topVolBranchSaddleIsoValue = ds.GetField("TopVolumeBranchSaddleIsoValue")
.GetData()
.AsArrayHandle<vtkm::cont::ArrayHandle<vtkm::Float32>>()
.ReadPortal();
vtkm::Id nSelectedBranches = topVolBranchGRId.GetNumberOfValues();
Edge expectedGRIdVolumeAtBranch0(38, 6);
Edge expectedEpsilonIsoAtBranch0(1, 50);
Edge expectedGRIdVolumeAtBranch1(50, 2);
Edge expectedEpsilonIsoAtBranch1(-1, 30);
for (vtkm::Id branch = 0; branch < nSelectedBranches; ++branch)
{
bool failed = false;
Edge computedGRIdVolume(topVolBranchGRId.Get(branch), topVolBranchVolume.Get(branch));
Edge computedEpsilonIso(topVolBranchSaddleEpsilon.Get(branch),
(vtkm::Id)topVolBranchSaddleIsoValue.Get(branch));
switch (branch)
{
case 0:
failed = !(computedGRIdVolume == expectedGRIdVolumeAtBranch0);
failed = (failed || !(computedEpsilonIso == expectedEpsilonIsoAtBranch0));
break;
case 1:
failed = !(computedGRIdVolume == expectedGRIdVolumeAtBranch1);
failed = (failed || !(computedEpsilonIso == expectedEpsilonIsoAtBranch1));
break;
default:
VTKM_TEST_ASSERT(false);
}
if (failed)
{
std::vector<Edge> expectedGRIdVolume{ expectedGRIdVolumeAtBranch0,
expectedGRIdVolumeAtBranch1 };
std::vector<Edge> expectedEpsilonIso{ expectedEpsilonIsoAtBranch0,
expectedEpsilonIsoAtBranch1 };
std::cout << "Top Branch Volume Results:" << std::endl;
std::cout << "Computed Top Branch Volume:branch=" << branch << std::endl;
std::cout << computedGRIdVolume.low << std::setw(14) << computedGRIdVolume.high
<< std::setw(5) << computedEpsilonIso.low << std::setw(14)
<< computedEpsilonIso.high << std::endl;
std::cout << "Expected Top Branch Volume:branch=" << branch << std::endl;
std::cout << expectedGRIdVolume[branch].low << std::setw(14)
<< expectedGRIdVolume[branch].high << std::setw(5)
<< expectedEpsilonIso[branch].low << std::setw(14)
<< expectedEpsilonIso[branch].high << std::endl;
VTKM_TEST_FAIL("Top Branch Volume Computation Failed!");
}
}
}
std::cout << "Expected Contour Tree" << std::endl;
for (std::size_t i = 0; i < expected.size(); i++)
{
std::cout << std::setw(12) << expected[i].low << std::setw(14) << expected[i].high
<< std::endl;
}
VTKM_TEST_FAIL("Branch Decomposition Failed!");
std::cout << "Top Branch Volume: Results Match!" << std::endl;
}
}

@ -16,6 +16,7 @@ set(headers
vtkm_declare_headers(${headers})
add_subdirectory(branch_decomposition)
add_subdirectory(select_top_volume_contours)
add_subdirectory(contourtree)
add_subdirectory(contourtree_augmented)
add_subdirectory(contourtree_distributed)

@ -600,7 +600,7 @@ inline void HierarchicalVolumetricBranchDecomposer::CollectEndsOfBranches(
{
IdArrayType superarcGRId;
vtkm::cont::make_ArrayHandlePermutation(supernodes, globalRegularIds);
vtkm::worklet::contourtree_augmented::PermuteArray<vtkm::Id>(
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
globalRegularIds, supernodes, superarcGRId);
std::stringstream resultStream;
@ -615,7 +615,7 @@ inline void HierarchicalVolumetricBranchDecomposer::CollectEndsOfBranches(
using InArrayHandleType = std::decay_t<decltype(inArray)>;
using ValueType = typename InArrayHandleType::ValueType;
vtkm::cont::ArrayHandle<ValueType> superarcValue;
vtkm::worklet::contourtree_augmented::PermuteArray<ValueType>(
vtkm::worklet::contourtree_augmented::PermuteArrayWithRawIndex<InArrayHandleType>(
inArray, supernodes, superarcValue);
vtkm::worklet::contourtree_augmented::PrintValues(
"Data Values", superarcValue, -1, resultStream);
@ -627,7 +627,6 @@ inline void HierarchicalVolumetricBranchDecomposer::CollectEndsOfBranches(
"Intrinsic Volumes", intrinsicVolumes, -1, resultStream);
vtkm::worklet::contourtree_augmented::PrintIndices(
"Dependent Volumes", dependentVolumes, -1, resultStream);
resultStream << std::endl;
VTKM_LOG_S(vtkm::cont::LogLevel::Info, resultStream.str());
}
@ -690,27 +689,27 @@ inline void HierarchicalVolumetricBranchDecomposer::CollectEndsOfBranches(
// actualOuterNodeValues[i] = dataValues[supernodes[outerNodes[actualSuperarcs[i]]]];
// actualOuterNodeRegularIds[i] = globalRegularIds[supernodes[outerNodes[actualSuperarcs[i]]]];
// }
// Solution: PermuteArray helps allocate the space so no need for explicit allocation
// PermuteArray also calls MaskedIndex
// Solution: PermuteArrayWithMaskedIndex helps allocate the space so no need for explicit allocation
// PermuteArrayWithMaskedIndex also calls MaskedIndex
// IdArrayType, size: nActualSuperarcs
IdArrayType actualBranchRoots;
vtkm::worklet::contourtree_augmented::PermuteArray<vtkm::Id>(
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
branchRoots, actualSuperarcs, actualBranchRoots);
// IdArrayType, size: nActualSuperarcs
IdArrayType actualOuterNodes;
vtkm::worklet::contourtree_augmented::PermuteArray<vtkm::Id>(
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
outerNodes, actualSuperarcs, actualOuterNodes);
// IdArrayType, size: nActualSuperarcs
IdArrayType actualOuterNodeLocalIds;
vtkm::worklet::contourtree_augmented::PermuteArray<vtkm::Id>(
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
supernodes, actualOuterNodes, actualOuterNodeLocalIds);
// IdArrayType, size: nActualSuperarcs
IdArrayType actualOuterNodeRegularIds;
vtkm::worklet::contourtree_augmented::PermuteArray<vtkm::Id>(
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
globalRegularIds, actualOuterNodeLocalIds, actualOuterNodeRegularIds);
auto resolveArray = [&](const auto& inArray) {
@ -744,7 +743,7 @@ inline void HierarchicalVolumetricBranchDecomposer::CollectEndsOfBranches(
// This is the real superarc local ID after permutation
IdArrayType permutedActualSuperarcs;
vtkm::worklet::contourtree_augmented::PermuteArray<vtkm::Id>(
vtkm::worklet::contourtree_augmented::PermuteArrayWithMaskedIndex<vtkm::Id>(
actualSuperarcs, sortedSuperarcs, permutedActualSuperarcs);
#ifdef DEBUG_HIERARCHICAL_VOLUMETRIC_BRANCH_DECOMPOSER

@ -942,13 +942,13 @@ inline void ActiveGraph::DebugPrint(const char* message, const char* fileName, l
// Active Vertex Arrays
IdArrayType activeIndices;
PermuteArray<vtkm::Id>(this->GlobalIndex, this->ActiveVertices, activeIndices);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->GlobalIndex, this->ActiveVertices, activeIndices);
IdArrayType activeFirst;
PermuteArray<vtkm::Id>(this->FirstEdge, this->ActiveVertices, activeFirst);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->FirstEdge, this->ActiveVertices, activeFirst);
IdArrayType activeOutdegree;
PermuteArray<vtkm::Id>(this->Outdegree, this->ActiveVertices, activeOutdegree);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->Outdegree, this->ActiveVertices, activeOutdegree);
IdArrayType activeHyperarcs;
PermuteArray<vtkm::Id>(this->Hyperarcs, this->ActiveVertices, activeHyperarcs);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->Hyperarcs, this->ActiveVertices, activeHyperarcs);
std::cout << "Active Vertex Arrays - Size: " << this->ActiveVertices.GetNumberOfValues()
<< std::endl;
PrintHeader(this->ActiveVertices.GetNumberOfValues());
@ -961,9 +961,9 @@ inline void ActiveGraph::DebugPrint(const char* message, const char* fileName, l
// Full Edge Arrays
IdArrayType farIndices;
PermuteArray<vtkm::Id>(this->GlobalIndex, this->EdgeFar, farIndices);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->GlobalIndex, this->EdgeFar, farIndices);
IdArrayType nearIndices;
PermuteArray<vtkm::Id>(this->GlobalIndex, this->EdgeNear, nearIndices);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->GlobalIndex, this->EdgeNear, nearIndices);
std::cout << "Full Edge Arrays - Size: " << this->EdgeNear.GetNumberOfValues() << std::endl;
PrintHeader(this->EdgeFar.GetNumberOfValues());
PrintIndices("Near", this->EdgeNear);
@ -974,9 +974,9 @@ inline void ActiveGraph::DebugPrint(const char* message, const char* fileName, l
// Active Edge Arrays
IdArrayType activeFarIndices;
PermuteArray<vtkm::Id>(this->EdgeFar, this->ActiveEdges, activeFarIndices);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->EdgeFar, this->ActiveEdges, activeFarIndices);
IdArrayType activeNearIndices;
PermuteArray<vtkm::Id>(this->EdgeNear, this->ActiveEdges, activeNearIndices);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->EdgeNear, this->ActiveEdges, activeNearIndices);
std::cout << "Active Edge Arrays - Size: " << this->ActiveEdges.GetNumberOfValues()
<< std::endl;
PrintHeader(this->ActiveEdges.GetNumberOfValues());
@ -987,9 +987,9 @@ inline void ActiveGraph::DebugPrint(const char* message, const char* fileName, l
// Edge Sorter Array
IdArrayType sortedFarIndices;
PermuteArray<vtkm::Id>(this->EdgeFar, this->EdgeSorter, sortedFarIndices);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->EdgeFar, this->EdgeSorter, sortedFarIndices);
IdArrayType sortedNearIndices;
PermuteArray<vtkm::Id>(this->EdgeNear, this->EdgeSorter, sortedNearIndices);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->EdgeNear, this->EdgeSorter, sortedNearIndices);
std::cout << "Edge Sorter - Size: " << this->EdgeSorter.GetNumberOfValues() << std::endl;
PrintHeader(this->EdgeSorter.GetNumberOfValues());
PrintIndices("Edge Sorter", this->EdgeSorter);

@ -75,7 +75,9 @@ namespace contourtree_augmented
// permute routines
template <typename ValueType, typename ArrayType>
inline void PermuteArray(const ArrayType& input, IdArrayType& permute, ArrayType& output)
inline void PermuteArrayWithMaskedIndex(const ArrayType& input,
IdArrayType& permute,
ArrayType& output)
{ // permuteValues()
using transform_type =
vtkm::cont::ArrayHandleTransform<IdArrayType, MaskedIndexFunctor<ValueType>>;
@ -110,6 +112,38 @@ inline void PermuteArray(const ArrayType& input, IdArrayType& permute, ArrayType
} // permuteValues()
// permute value type arrays
template <typename ArrayType>
inline void PermuteArrayWithRawIndex(const ArrayType& input,
IdArrayType& permute,
ArrayType& output)
{ // PermuteArrayWithRawIndex()
using permute_type = vtkm::cont::ArrayHandlePermutation<IdArrayType, ArrayType>;
// resize the output, i.e., output.resize(permute.size());
vtkm::Id permNumValues = permute.GetNumberOfValues();
vtkm::Id outNumValues = output.GetNumberOfValues();
if (permNumValues > outNumValues)
{
output.Allocate(permNumValues);
}
else if (permNumValues < outNumValues)
{
output.Allocate(permNumValues, vtkm::CopyFlag::On);
} // else the output has already the correct size
// The following is equivalent to doing the following in serial
//
// for (vtkm::Id entry = 0; entry < permute.size(); entry++)
// output[entry] = input[permute[entry]];
//
// fancy vtkm array so that we do not actually copy any data here
permute_type permutedInput(permute, input);
// Finally, copy the permuted values to the output array
vtkm::cont::ArrayCopyDevice(permutedInput, output);
} // PermuteArrayWithRawIndex()
// transform functor used in ContourTreeMesh to flag indicies as other when using the CombinedVectorClass
struct MarkOther
{

@ -268,13 +268,13 @@ inline void ContourTreeMaker::ComputeHyperAndSuperStructure()
// we have to permute a bunch of arrays, so let's have some temporaries to store them
IdArrayType permutedHyperparents;
PermuteArray<vtkm::Id>(
PermuteArrayWithMaskedIndex<vtkm::Id>(
this->ContourTreeResult.Hyperparents, this->ContourTreeResult.Hypernodes, permutedHyperparents);
IdArrayType permutedSupernodes;
PermuteArray<vtkm::Id>(
PermuteArrayWithMaskedIndex<vtkm::Id>(
this->ContourTreeResult.Supernodes, this->ContourTreeResult.Hypernodes, permutedSupernodes);
IdArrayType permutedSuperarcs;
PermuteArray<vtkm::Id>(
PermuteArrayWithMaskedIndex<vtkm::Id>(
this->ContourTreeResult.Superarcs, this->ContourTreeResult.Hypernodes, permutedSuperarcs);
// now we establish the reverse index array
@ -284,7 +284,7 @@ inline void ContourTreeMaker::ComputeHyperAndSuperStructure()
// for (vtkm::Id supernode = 0; supernode < this->ContourTreeResult.Supernodes.size(); supernode++)
// superSortIndex[this->ContourTreeResult.Hypernodes[supernode]] = supernode;
//typedef vtkm::cont::ArrayHandlePermutation<IdArrayType, vtkm::cont::ArrayHandleIndex> PermuteArrayHandleIndex;
//typedef vtkm::cont::ArrayHandlePermutation<IdArrayType, vtkm::cont::ArrayHandleIndex> PermuteArrayWithMaskedIndexHandleIndex;
vtkm::cont::ArrayHandlePermutation<IdArrayType, IdArrayType> permutedSuperSortIndex(
this->ContourTreeResult.Hypernodes, // index array
superSortIndex); // value array
@ -307,7 +307,7 @@ inline void ContourTreeMaker::ComputeHyperAndSuperStructure()
// we will permute the hyperarcs & copy them back with the new supernode target IDs
IdArrayType permutedHyperarcs;
PermuteArray<vtkm::Id>(
PermuteArrayWithMaskedIndex<vtkm::Id>(
this->ContourTreeResult.Hyperarcs, this->ContourTreeResult.Hypernodes, permutedHyperarcs);
contourtree_maker_inc_ns::ComputeHyperAndSuperStructure_PermuteArcs permuteHyperarcsWorklet;
this->Invoke(
@ -315,9 +315,9 @@ inline void ContourTreeMaker::ComputeHyperAndSuperStructure()
// now swizzle the WhenTransferred value
IdArrayType permutedWhenTransferred;
PermuteArray<vtkm::Id>(this->ContourTreeResult.WhenTransferred,
this->ContourTreeResult.Hypernodes,
permutedWhenTransferred);
PermuteArrayWithMaskedIndex<vtkm::Id>(this->ContourTreeResult.WhenTransferred,
this->ContourTreeResult.Hypernodes,
permutedWhenTransferred);
vtkm::cont::Algorithm::Copy(permutedWhenTransferred, this->ContourTreeResult.WhenTransferred);
// now we compress both the hypernodes & Hyperarcs

@ -72,7 +72,7 @@ namespace contourtree_augmented
{
// local constants to allow changing the spacing as needed
constexpr int PRINT_WIDTH = 12;
constexpr int PRINT_WIDTH = 14;
constexpr int PREFIX_WIDTH = 30;
template <typename T, typename StorageType>

@ -659,11 +659,12 @@ public:
IdArrayType permutedBranches;
permutedBranches.Allocate(nSupernodes);
PermuteArray<vtkm::Id>(whichBranch, supernodeSorter, permutedBranches);
PermuteArrayWithMaskedIndex<vtkm::Id>(whichBranch, supernodeSorter, permutedBranches);
IdArrayType permutedRegularID;
permutedRegularID.Allocate(nSupernodes);
PermuteArray<vtkm::Id>(contourTree.Supernodes, supernodeSorter, permutedRegularID);
PermuteArrayWithMaskedIndex<vtkm::Id>(
contourTree.Supernodes, supernodeSorter, permutedRegularID);
#ifdef DEBUG_PRINT
std::cout << "VI A. Sorted into Branches" << std::endl;

@ -0,0 +1,146 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#ifndef vtk_m_filter_scalar_topology_worklet_branch_decomposition_select_top_volume_contours_BranchVolumeComparator_h
#define vtk_m_filter_scalar_topology_worklet_branch_decomposition_select_top_volume_contours_BranchVolumeComparator_h
#include <vtkm/worklet/WorkletMapField.h>
namespace vtkm
{
namespace worklet
{
namespace scalar_topology
{
namespace select_top_volume_contours
{
using IdArrayType = vtkm::worklet::contourtree_augmented::IdArrayType;
// Implementation of BranchVolumeComparator
class BranchVolumeComparatorImpl
{
public:
using IdPortalType = typename IdArrayType::ReadPortalType;
// constructor
VTKM_CONT
BranchVolumeComparatorImpl(const IdArrayType& branchRoots,
const IdArrayType& branchVolume,
vtkm::cont::DeviceAdapterId device,
vtkm::cont::Token& token)
: branchRootsPortal(branchRoots.PrepareForInput(device, token))
, branchVolumePortal(branchVolume.PrepareForInput(device, token))
{ // constructor
} // constructor
// () operator - gets called to do comparison
VTKM_EXEC
bool operator()(const vtkm::Id& i, const vtkm::Id& j) const
{ // operator()
vtkm::Id volumeI = this->branchVolumePortal.Get(i);
vtkm::Id volumeJ = this->branchVolumePortal.Get(j);
// primary sort on branch volume
if (volumeI > volumeJ)
return true;
if (volumeI < volumeJ)
return false;
vtkm::Id branchI =
vtkm::worklet::contourtree_augmented::MaskedIndex(this->branchRootsPortal.Get(i));
vtkm::Id branchJ =
vtkm::worklet::contourtree_augmented::MaskedIndex(this->branchRootsPortal.Get(j));
// secondary sort on branch ID
return (branchI < branchJ);
} // operator()
private:
IdPortalType branchRootsPortal;
IdPortalType branchVolumePortal;
}; // BranchVolumeComparatorImpl
/// <summary>
/// Comparator of branch volume. Higher volume comes first
/// </summary>
class BranchVolumeComparator : public vtkm::cont::ExecutionObjectBase
{
public:
// constructor
VTKM_CONT
BranchVolumeComparator(const IdArrayType& branchRoots, const IdArrayType& branchVolume)
: BranchRoots(branchRoots)
, BranchVolume(branchVolume)
{
}
VTKM_CONT BranchVolumeComparatorImpl PrepareForExecution(vtkm::cont::DeviceAdapterId device,
vtkm::cont::Token& token) const
{
return BranchVolumeComparatorImpl(this->BranchRoots, this->BranchVolume, device, token);
}
private:
IdArrayType BranchRoots;
IdArrayType BranchVolume;
}; // BranchVolumeComparator
} // namespace select_top_volume_contours
} // namespace scalar_topology
} // namespace worklet
} // namespace vtkm
#endif

@ -0,0 +1,19 @@
##============================================================================
## 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.
##============================================================================
set(headers
ClarifyBranchEndSupernodeTypeWorklet.h
UpdateInfoByBranchDirectionWorklet.h
GetBranchVolumeWorklet.h
BranchVolumeComparator.h
)
#-----------------------------------------------------------------------------
vtkm_declare_headers(${headers})

@ -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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#ifndef vtk_m_filter_scalar_topology_worklet_select_top_volume_contours_ClarifyBranchEndSupernodeTypeWorklet_h
#define vtk_m_filter_scalar_topology_worklet_select_top_volume_contours_ClarifyBranchEndSupernodeTypeWorklet_h
#include <vtkm/worklet/WorkletMapField.h>
namespace vtkm
{
namespace worklet
{
namespace scalar_topology
{
namespace select_top_volume_contours
{
// For special branches that only have one superarc
// Clarify which end is leaf and which is saddle
class ClarifyBranchEndSupernodeTypeWorklet : public vtkm::worklet::WorkletMapField
{
public:
using ControlSignature = void(
FieldIn lowerSuperarcId, // (input) lower end superarc ID
FieldIn lowerIntrinsic, // (input) lower end superarc intrisic volume
FieldIn upperSuperarcId, // (input) upper end superarc ID
FieldIn upperIntrinsic, // (input) upper end superarc intrisic volume
FieldIn branchRoot, // (input) branch root superarc ID
FieldInOut isLowerLeaf, // (input/output) bool, whether the lower end is a leaf
FieldInOut isUpperLeaf // (input/output) bool, whether the upper end is a leaf
);
using ExecutionSignature = void(_1, _2, _3, _4, _5, _6, _7);
using InputDomain = _1;
/// Constructor
VTKM_EXEC_CONT
ClarifyBranchEndSupernodeTypeWorklet(const vtkm::Id tVol)
: totalVolume(tVol)
{
}
/// The functor checks the direction of the branch
VTKM_EXEC void operator()(const vtkm::Id& lowerSuperarcId,
const vtkm::Id& lowerIntrinsic,
const vtkm::Id& upperSuperarcId,
const vtkm::Id& upperIntrinsic,
const vtkm::Id& branchRoot,
bool& isLowerLeaf,
bool& isUpperLeaf) const
{
// do nothing: not a "leaf-leaf" branch
if (!isLowerLeaf || !isUpperLeaf)
return;
// do nothing: actual leaf-leaf branch
if (lowerIntrinsic == totalVolume - 1 && lowerIntrinsic == upperIntrinsic)
return;
// do something: fake leaf-leaf branch
// we already exclude the case of only one superarc above
vtkm::Id maskedLowerId = vtkm::worklet::contourtree_augmented::MaskedIndex(lowerSuperarcId);
vtkm::Id maskedUpperId = vtkm::worklet::contourtree_augmented::MaskedIndex(upperSuperarcId);
vtkm::Id maskedRoot = vtkm::worklet::contourtree_augmented::MaskedIndex(branchRoot);
if (maskedLowerId == maskedRoot && maskedUpperId == maskedRoot)
{
const bool isAscending = vtkm::worklet::contourtree_augmented::IsAscending(lowerSuperarcId);
if (isAscending)
isUpperLeaf = false;
else
isLowerLeaf = false;
}
}
private:
const vtkm::Id totalVolume;
}; // ClarifyBranchEndSupernodeTypeWorklet
} // namespace select_top_volume_contours
} // namespace scalar_topology
} // namespace worklet
} // namespace vtkm
#endif

@ -0,0 +1,137 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#ifndef vtk_m_filter_scalar_topology_worklet_select_top_volume_contours_GetBranchVolumeWorklet_h
#define vtk_m_filter_scalar_topology_worklet_select_top_volume_contours_GetBranchVolumeWorklet_h
#include <vtkm/worklet/WorkletMapField.h>
namespace vtkm
{
namespace worklet
{
namespace scalar_topology
{
namespace select_top_volume_contours
{
using IdArrayType = vtkm::worklet::contourtree_augmented::IdArrayType;
/// <summary>
/// worklet to check the direction of branch
/// return true if the branch inner superarc points to the senior-most node
/// return false if the branch inner superarc comes from the senior-most node
/// trick: if the branch is the main branch, we return true for computation later
/// </summary>
class GetBranchVolumeWorklet : public vtkm::worklet::WorkletMapField
{
public:
using ControlSignature = void(
FieldIn lowerDirection, // (input) lower end superarc ID with direction information
FieldIn lowerIntrinsic, // (input) lower end superarc intrisic volume
FieldIn lowerDependent, // (input) lower end superarc dependent volume
FieldIn upperDirection, // (input) upper end superarc ID with direction information
FieldIn upperIntrinsic, // (input) upper end superarc intrisic volume
FieldIn upperDependent, // (input) upper end superarc dependent volume
FieldIn isLowerLeaf, // (input) bool, whether the lower end is a leaf
FieldIn isUpperLeaf, // (input) bool, whether the upper end is a leaf
FieldOut branchVolume // (output) volume of the branch
);
using ExecutionSignature = _9(_1, _2, _3, _4, _5, _6, _7, _8);
using InputDomain = _1;
/// Constructor
VTKM_EXEC_CONT
GetBranchVolumeWorklet(const vtkm::Id tVol)
: totalVolume(tVol)
{
}
/// The functor checks the direction of the branch
VTKM_EXEC vtkm::Id operator()(const vtkm::Id& lowerDirection,
const vtkm::Id& lowerIntrinsic,
const vtkm::Id& lowerDependent,
const vtkm::Id& upperDirection,
const vtkm::Id& upperIntrinsic,
const vtkm::Id& upperDependent,
const bool& isLowerLeaf,
const bool& isUpperLeaf) const
{
if (isLowerLeaf && isUpperLeaf)
return totalVolume;
// if the branch is a minimum-saddle branch
// if the upper end superarc direction is pointing up, then dependent; otherwise, reverse
if (isLowerLeaf)
return contourtree_augmented::IsAscending(upperDirection)
? upperDependent
: totalVolume - upperDependent + upperIntrinsic;
// if the branch is a maximum-saddle branch
// if the lower end superarc direction is pointing down, then true; otherwise, false
if (isUpperLeaf)
return !contourtree_augmented::IsAscending(lowerDirection)
? lowerDependent
: totalVolume - lowerDependent + lowerIntrinsic;
// in case of fallout, should never reach
return 0;
}
private:
const vtkm::Id totalVolume;
}; // GetBranchVolumeWorklet
} // namespace select_top_volume_contours
} // namespace scalar_topology
} // namespace worklet
} // namespace vtkm
#endif

@ -0,0 +1,124 @@
//============================================================================
// 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 (c) 2018, The Regents of the University of California, through
// Lawrence Berkeley National Laboratory (subject to receipt of any required approvals
// from the U.S. Dept. of Energy). All rights reserved.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// (1) Redistributions of source code must retain the above copyright notice, this
// list of conditions and the following disclaimer.
//
// (2) Redistributions in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// (3) Neither the name of the University of California, Lawrence Berkeley National
// Laboratory, U.S. Dept. of Energy nor the names of its contributors may be
// used to endorse or promote products derived from this software without
// specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
// IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
// INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
// BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
// DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
// OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED
// OF THE POSSIBILITY OF SUCH DAMAGE.
//
//=============================================================================
//
// This code is an extension of the algorithm presented in the paper:
// Parallel Peak Pruning for Scalable SMP Contour Tree Computation.
// Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens.
// Proceedings of the IEEE Symposium on Large Data Analysis and Visualization
// (LDAV), October 2016, Baltimore, Maryland.
//
// The PPP2 algorithm and software were jointly developed by
// Hamish Carr (University of Leeds), Gunther H. Weber (LBNL), and
// Oliver Ruebel (LBNL)
//==============================================================================
#ifndef vtk_m_filter_scalar_topology_worklet_select_top_volume_contours_UpdateInfoByBranchDirectionWorklet_h
#define vtk_m_filter_scalar_topology_worklet_select_top_volume_contours_UpdateInfoByBranchDirectionWorklet_h
#include <vtkm/worklet/WorkletMapField.h>
namespace vtkm
{
namespace worklet
{
namespace scalar_topology
{
namespace select_top_volume_contours
{
/// <summary>
/// 1. Get the saddle end isovalue
/// 2. Get Epsilon direction near the branch saddle end
/// If main branch, epsilon is 0.
/// Otherwise, -1 if the branch is the lower leaf branch, or 1 if upper leaf branch.
/// </summary>
template <typename ValueType>
class UpdateInfoByBranchDirectionWorklet : public vtkm::worklet::WorkletMapField
{
public:
using ControlSignature = void(
FieldIn isLowerLeaf, // (input) bool, whether the lower end is a leaf
FieldIn isUpperLeaf, // (input) bool, whether the upper end is a leaf
FieldIn upperEndValue, // (input) value type, the branch upper end value
FieldIn lowerEndValue, // (input) value type, the branch lower end value
FieldOut saddleEpsilon, // (output) vtkm::Id, epsilon direction around the saddle isovalue
FieldOut saddleValue // (output) value type, the saddle isovalue
);
using ExecutionSignature = void(_1, _2, _3, _4, _5, _6);
using InputDomain = _1;
/// Constructor
VTKM_EXEC_CONT
UpdateInfoByBranchDirectionWorklet()
{ // constructor
} // constructor
/// The functor returns the isovalue and the epsilon direction around the saddle end of the branch
VTKM_EXEC void operator()(const bool isLowerLeaf,
const bool isUpperLeaf,
const ValueType upperEndValue,
const ValueType lowerEndValue,
vtkm::Id& saddleEpsilon,
ValueType& saddleValue) const
{
// NOTE: for the main branch, the saddle value is undefined,
// because both upper and lower ends are leaf nodes.
// Let's use upperEndValue here to make the output not random.
if (isLowerLeaf && isUpperLeaf)
{
saddleEpsilon = 0;
saddleValue = upperEndValue;
}
else
{
saddleEpsilon = isLowerLeaf ? -1 : 1;
saddleValue = isLowerLeaf ? upperEndValue : lowerEndValue;
}
}
}; // EpsilonFromBranchDirection
} // namespace select_top_volume_contours
} // namespace scalar_topology
} // namespace worklet
} // namespace vtkm
#endif