Commit bc550a57 authored by Cory Quammen's avatar Cory Quammen Committed by Kitware Robot
Browse files

Merge topic '18514-regionid-ordering-options'

0ca613af Fix compiler warnings
0502cdf4 Add CELL_COUNT_ASCENDING RegionIdAssignmentMode
3ae247c9 Update ParallelConnectivity test to exercise RegionIdAssignmentMode
be24b724 Add support for reordering RegionIds by RegionIdAssignmentMode
97ba574f

 Add RegionId assignment mode

Acked-by: default avatarKitware Robot <kwrobot@kitware.com>
Acked-by: default avatarUtkarsh Ayachit <utkarsh.ayachit@kitware.com>
Merge-request: !4812
parents 54a89e66 0ca613af
......@@ -31,6 +31,8 @@
#include "vtkUnstructuredGrid.h"
#include "vtkIdTypeArray.h"
#include <map>
vtkObjectFactoryNewMacro(vtkConnectivityFilter);
// Construct with default extraction mode to extract largest regions.
......@@ -39,6 +41,7 @@ vtkConnectivityFilter::vtkConnectivityFilter()
this->RegionSizes = vtkIdTypeArray::New();
this->ExtractionMode = VTK_EXTRACT_LARGEST_REGION;
this->ColorRegions = 0;
this->RegionIdAssignmentMode = UNSPECIFIED;
this->ScalarConnectivity = 0;
this->ScalarRange[0] = 0.0;
......@@ -350,6 +353,8 @@ int vtkConnectivityFilter::RequestData(
// if coloring regions; send down new scalar data
if ( this->ColorRegions )
{
this->OrderRegionIds(this->NewScalars, this->NewCellScalars);
int idx = outputPD->AddArray(this->NewScalars);
outputPD->SetActiveAttribute(idx, vtkDataSetAttributes::SCALARS);
idx = outputCD->AddArray(this->NewCellScalars);
......@@ -613,6 +618,75 @@ void vtkConnectivityFilter::TraverseAndMark (vtkDataSet *input)
} //while wave is not empty
}
void vtkConnectivityFilter::OrderRegionIds(vtkIdTypeArray* pointRegionIds, vtkIdTypeArray* cellRegionIds)
{
if (this->ColorRegions)
{
if (this->RegionIdAssignmentMode == CELL_COUNT_DESCENDING ||
this->RegionIdAssignmentMode == CELL_COUNT_ASCENDING)
{
// Use a multimap to handle cases where more than one region has the same number of cells.
std::multimap<vtkIdType, vtkIdType> cellCountToRegionId;
typedef std::multimap<vtkIdType, vtkIdType>::value_type ValueType;
vtkIdType numRegions = this->RegionSizes->GetNumberOfTuples();
for (vtkIdType regionId = 0; regionId < numRegions; ++regionId)
{
ValueType value(this->RegionSizes->GetValue(regionId), regionId);
cellCountToRegionId.insert(value);
}
// Now reverse iterate through the sorted multimap to process the RegionIds
// from largest to smallest and create a map from the old RegionId to the new
// RegionId
std::map<vtkIdType, vtkIdType> oldToNew;
vtkIdType counter = 0;
if (this->RegionIdAssignmentMode == CELL_COUNT_ASCENDING)
{
for (auto iter = cellCountToRegionId.begin(); iter != cellCountToRegionId.end(); ++iter)
{
auto regionCount = iter->first;
auto regionId = iter->second;
// Re-order the region sizes based on the sorting
this->RegionSizes->SetValue(counter, regionCount);
// Create map from old to new RegionId
oldToNew[regionId] = counter++;
}
}
else // CELL_COUNT_DESCENDING
{
for (auto iter = cellCountToRegionId.rbegin(); iter != cellCountToRegionId.rend(); ++iter)
{
auto regionCount = iter->first;
auto regionId = iter->second;
// Re-order the region sizes based on the sorting
this->RegionSizes->SetValue(counter, regionCount);
// Create map from old to new RegionId
oldToNew[regionId] = counter++;
}
}
vtkIdType numPts = pointRegionIds->GetNumberOfTuples();
for (vtkIdType i = 0; i < numPts; ++i)
{
vtkIdType oldValue = pointRegionIds->GetValue(i);
pointRegionIds->SetValue(i, oldToNew[oldValue]);
}
vtkIdType numCells = cellRegionIds->GetNumberOfTuples();
for (vtkIdType i = 0; i < numCells; ++i)
{
vtkIdType oldValue = cellRegionIds->GetValue(i);
cellRegionIds->SetValue(i, oldToNew[oldValue]);
}
}
// else UNSPECIFIED mode
}
}
// Obtain the number of connected regions.
int vtkConnectivityFilter::GetNumberOfExtractedRegions()
{
......
......@@ -180,6 +180,24 @@ public:
vtkBooleanMacro(ColorRegions,vtkTypeBool);
//@}
/**
* Enumeration of the various ways to assign RegionIds when
* the ColorRegions option is on.
*/
enum RegionIdAssignment {
UNSPECIFIED,
CELL_COUNT_DESCENDING,
CELL_COUNT_ASCENDING
};
//@{
/**
* Set/get mode controlling how RegionIds are assigned.
*/
//@}
vtkSetMacro(RegionIdAssignmentMode, int);
vtkGetMacro(RegionIdAssignmentMode, int);
//@{
/**
* Set/get the desired precision for the output types. See the documentation
......@@ -215,8 +233,12 @@ protected:
vtkTypeBool ScalarConnectivity;
double ScalarRange[2];
int RegionIdAssignmentMode;
void TraverseAndMark(vtkDataSet *input);
void OrderRegionIds(vtkIdTypeArray* pointRegionIds, vtkIdTypeArray* cellRegionIds);
private:
// used to support algorithm execution
vtkFloatArray *CellScalars;
......
......@@ -15,9 +15,11 @@
#include "vtkConnectivityFilter.h"
#include "vtkCellData.h"
#include "vtkContourFilter.h"
#include "vtkDataSetTriangleFilter.h"
#include "vtkDistributedDataFilter.h"
#include "vtkIdTypeArray.h"
#include "vtkMPIController.h"
#include "vtkPUnstructuredGridGhostCellsGenerator.h"
#include "vtkPConnectivityFilter.h"
......@@ -98,13 +100,97 @@ int RunParallelConnectivity(const char* fname, vtkAlgorithm::DesiredOutputPrecis
returnValue = EXIT_FAILURE;
}
// Check that assigning RegionIds by number of cells (descending) works
connectivity->SetRegionIdAssignmentMode(vtkConnectivityFilter::CELL_COUNT_DESCENDING);
connectivity->ColorRegionsOn();
connectivity->SetExtractionModeToAllRegions();
removeGhosts->Update();
numberOfRegions = connectivity->GetNumberOfExtractedRegions();
vtkPointSet* ghostOutput = vtkPointSet::SafeDownCast(removeGhosts->GetOutput());
vtkIdType numberOfCells = ghostOutput->GetNumberOfCells();
vtkIdType globalNumberOfCells = 0;
contr->AllReduce(&numberOfCells, &globalNumberOfCells, 1, vtkCommunicator::SUM_OP);
std::vector<vtkIdType> regionCounts(connectivity->GetNumberOfExtractedRegions(), 0);
// Count up cells with RegionIds
auto regionIdArray = vtkIdTypeArray::SafeDownCast(ghostOutput->GetCellData()->GetArray("RegionId"));
for (vtkIdType cellId = 0; cellId < numberOfCells; ++cellId)
{
vtkIdType regionId = regionIdArray->GetValue(cellId);
regionCounts[regionId]++;
}
// Sum up region counts across processes
std::vector<vtkIdType> globalRegionCounts(regionCounts.size(), 0);
contr->AllReduce(regionCounts.data(), globalRegionCounts.data(),
static_cast<vtkIdType>(regionCounts.size()), vtkCommunicator::SUM_OP);
if (me == 0)
{
bool printCounts = false;
for (vtkIdType i = 1; i < numberOfRegions; ++i)
{
if (globalRegionCounts[i] > globalRegionCounts[i-1])
{
std::cerr << "Region " << i << " is larger than region " << i-1 << std::endl;
printCounts = true;
returnValue = EXIT_FAILURE;
break;
}
}
if (printCounts)
{
for (vtkIdType i = 0; i < numberOfRegions; ++i)
{
std::cout << "Region " << i << " has " << globalRegionCounts[i] << " cells" << std::endl;
}
}
}
// Check that assignment RegionIds by number of cells (ascending) works
connectivity->SetRegionIdAssignmentMode(vtkConnectivityFilter::CELL_COUNT_ASCENDING);
removeGhosts->Update();
std::fill(regionCounts.begin(), regionCounts.end(), 0);
regionIdArray = vtkIdTypeArray::SafeDownCast(ghostOutput->GetCellData()->GetArray("RegionId"));
for (vtkIdType cellId = 0; cellId < numberOfCells; ++cellId)
{
vtkIdType regionId = regionIdArray->GetValue(cellId);
regionCounts[regionId]++;
}
// Sum up region counts across processes
globalRegionCounts = std::vector<vtkIdType>(regionCounts.size(), 0);
contr->AllReduce(regionCounts.data(), globalRegionCounts.data(),
static_cast<vtkIdType>(regionCounts.size()), vtkCommunicator::SUM_OP);
if (me == 0)
{
bool printCounts = false;
for (vtkIdType i = 1; i < numberOfRegions; ++i)
{
if (globalRegionCounts[i] < globalRegionCounts[i-1])
{
std::cerr << "Region " << i << " is smaller than " << i-1 << std::endl;
printCounts = true;
returnValue = EXIT_FAILURE;
break;
}
}
if (printCounts)
{
for (vtkIdType i = 0; i < numberOfRegions; ++i)
{
std::cout << "Region " << i << " has " << globalRegionCounts[i] << " cells" << std::endl;
}
}
}
// Check the number of cells in the largest region when the extraction mode
// is set to largest region.
connectivity->SetExtractionModeToLargestRegion();
removeGhosts->Update();
int numberOfCells =
numberOfCells =
vtkPointSet::SafeDownCast(removeGhosts->GetOutput())->GetNumberOfCells();
int globalNumberOfCells = 0;
globalNumberOfCells = 0;
contr->AllReduce(&numberOfCells, &globalNumberOfCells, 1, vtkCommunicator::SUM_OP);
int expectedNumberOfCells = 2124;
......
......@@ -482,11 +482,13 @@ int vtkPConnectivityFilter::RequestData(
int saveScalarConnectivity = this->ScalarConnectivity;
int saveExtractionMode = this->ExtractionMode;
int saveColorRegions = this->ColorRegions;
int saveRegionIdAssignmentMode = this->RegionIdAssignmentMode;
// Overwrite custom member variables temporarily.
this->ScalarConnectivity = 0;
this->ExtractionMode = VTK_EXTRACT_ALL_REGIONS;
this->ColorRegions = 1;
this->RegionIdAssignmentMode = UNSPECIFIED;
// Invoke the connectivity algorithm in the superclass.
success = this->Superclass::RequestData(request, inputVector, outputVector);
......@@ -494,6 +496,7 @@ int vtkPConnectivityFilter::RequestData(
this->ScalarConnectivity = saveScalarConnectivity;
this->ExtractionMode = saveExtractionMode;
this->ColorRegions = saveColorRegions;
this->RegionIdAssignmentMode = saveRegionIdAssignmentMode;
}
else
{
......@@ -794,9 +797,15 @@ int vtkPConnectivityFilter::RequestData(
std::vector< vtkIdType > localRegionSizes(numContiguousLabels, 0);
if (cellRegionIds)
{
// Iterate over cells and count how many are in different regions.
// Iterate over cells and count how many are in different regions. Count only non-ghost cells.
vtkUnsignedCharArray* cellGhostArray = output->GetCellGhostArray();
for (vtkIdType i = 0; i < cellRegionIds->GetNumberOfValues(); ++i)
{
if (cellGhostArray && cellGhostArray->GetTypedComponent(i, 0) &
vtkDataSetAttributes::DUPLICATECELL)
{
continue;
}
localRegionSizes[cellRegionIds->GetValue(i)]++;
}
}
......@@ -815,15 +824,26 @@ int vtkPConnectivityFilter::RequestData(
this->RegionSizes->SetTypedTuple(i, &globalRegionSizes[i]);
}
// Potentially reorder RegionIds in the output arrays.
this->OrderRegionIds(pointRegionIds, cellRegionIds);
if (this->ExtractionMode == VTK_EXTRACT_LARGEST_REGION ||
this->ExtractionMode == VTK_EXTRACT_CLOSEST_POINT_REGION)
{
double threshold = 0.0;
if (this->ExtractionMode == VTK_EXTRACT_LARGEST_REGION)
{
vtkIdType largestRegionId =
std::distance(globalRegionSizes.begin(),
std::max_element(globalRegionSizes.begin(), globalRegionSizes.end()));
vtkIdType largestRegionCount = 0;
vtkIdType largestRegionId = 0;
for (vtkIdType i = 0; i < this->RegionSizes->GetNumberOfTuples(); ++i)
{
vtkIdType candidateCount = this->RegionSizes->GetValue(i);
if (candidateCount > largestRegionCount)
{
largestRegionCount = candidateCount;
largestRegionId = i;
}
}
threshold = largestRegionId;
}
else if (this->ExtractionMode == VTK_EXTRACT_CLOSEST_POINT_REGION)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment