Commit 10b1c2ef authored by Cory Quammen's avatar Cory Quammen

Add vtkDataObjectTreeToPointSet filter

If the ExtractSurfaces option is on, this filter produces a single
vtkPolyData output from the surfaces extracted from blocks in the
input. If off, the filter produces a single vtkUnstructuredGrid output
from the entire blocks in the input.

This class is now used for the "Merge Blocks" filter in ParaView.
parent 43d03e63
......@@ -13163,22 +13163,39 @@ Every Nth Point.
<!-- End of ConvertSelection -->
</SourceProxy>
<!-- ==================================================================== -->
<SourceProxy class="vtkCompositeDataToUnstructuredGridFilter"
<SourceProxy class="vtkDataObjectTreeToPointSetFilter"
label="Merge Blocks"
name="MergeBlocks">
<Documentation short_help="Appends vtkCompositeDataSet leaves into a single vtkUnstructuredGrid">
vtkCompositeDataToUnstructuredGridFilter appends all vtkDataSet leaves of
the input composite dataset to a single unstructured grid. The subtree to
be combined can be chosen using the SubTreeCompositeIndex. If the
SubTreeCompositeIndex is a leaf node, then no appending is
required.</Documentation>
<Documentation short_help="Appends vtkCompositeDataSet leaves into a single vtkUnstructuredGrid or vtkPolyData">
Merge Blocks appends all vtkDataSet leaves of the input composite dataset to a single
unstructured grid (or polydata if all leaves are polydata). The subtree to
be combined can be chosen using the SubTreeCompositeIndex. If the SubTreeCompositeIndex
is a leaf node, then no appending is required nor performed.</Documentation>
<InputProperty command="SetInputConnection"
name="Input">
<DataTypeDomain name="input_type">
<DataType value="vtkCompositeDataSet" />
<DataType value="vtkDataObjectTree" />
</DataTypeDomain>
<Documentation>Set the input composite dataset.</Documentation>
</InputProperty>
<IntVectorProperty command="SetOutputDataSetType"
default_values="4"
number_of_elements="1"
name="OutputDataSetType"
panel_visibility="advanced">
<EnumerationDomain name="enum">
<Entry text="Polygonal Mesh"
value="0" />
<Entry text="Unstructured Grid"
value="4" />
</EnumerationDomain>
<Documentation>Determines the output type produced by this filter. Only blocks compatible
with the output type will be merged in the output. For example, if the output type is
"Polygonal Mesh", then blocks of type "Image Data", "Structured Grid", etc. will not be merged.
If the output type is "UnstructuredGrid", then blocks of any type will be merged in the output.
Defaults to "Unstructured Grid".
</Documentation>
</IntVectorProperty>
<IntVectorProperty command="SetSubTreeCompositeIndex"
default_values="0"
is_internal="1"
......@@ -13200,14 +13217,34 @@ Every Nth Point.
name="MergePoints"
number_of_elements="1">
<BooleanDomain name="bool" />
<Documentation>Merge points within a distance specified by the **Tolerance**
property.</Documentation>
</IntVectorProperty>
<DoubleVectorProperty command="SetTolerance"
default_values="0"
label="Tolerance"
name="Tolerance"
number_of_elements="1"
panel_visibility="advanced">
<Documentation>Set the tolerance for merging points if **Merge Points** is enabled.</Documentation>
<Hints>
<PropertyWidgetDecorator type="GenericDecorator"
mode="visibility"
property="MergePoints"
value="1" />
<!-- show this widget when MergePoints==1 and showing advanced options -->
</Hints>
</DoubleVectorProperty>
<IntVectorProperty command="SetToleranceIsAbsolute"
default_values="0"
label="Tolerance"
name="Tolerance"
name="ToleranceIsAbsolute"
number_of_elements="1"
panel_visibility="advanced">
<Documentation>Set the tolerance for merging points if **Merge Points** is enabled.</Documentation>
<BooleanDomain name="bool" />
<Documentation>This property determines whether to treat the **Tolerance**
property as absolute (points closer than Tolerance are merged) or relative
(points closer than a fraction of the input data's bounding box diagonal length)
when performing point merging.</Documentation>
<Hints>
<PropertyWidgetDecorator type="GenericDecorator"
mode="visibility"
......@@ -13215,7 +13252,7 @@ Every Nth Point.
value="1" />
<!-- show this widget when MergePoints==1 and showing advanced options -->
</Hints>
</DoubleVectorProperty>
</IntVectorProperty>
<!-- End of MergeBlocks -->
</SourceProxy>
<!-- ==================================================================== -->
......
......@@ -39,6 +39,7 @@ set(classes
vtkCleanArrays
vtkCompositeDataToUnstructuredGridFilter
vtkContext2DScalarBarActor
vtkDataObjectTreeToPointSetFilter
vtkImageCompressor
vtkImageTransparencyFilter
vtkKdTreeGenerator
......
/*=========================================================================
Program: ParaView
Module: vtkDataObjectTreeToPointSetFilter.cxx
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.paraview.org/HTML/Copyright.html 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 "vtkDataObjectTreeToPointSetFilter.h"
#include "vtkAppendDataSets.h"
#include "vtkAppendFilter.h" // For unstructured grid output
#include "vtkAppendPolyData.h" // For polydata output
#include "vtkCleanArrays.h"
#include "vtkDataObjectTree.h"
#include "vtkDataObjectTreeIterator.h"
#include "vtkDataObjectTypes.h"
#include "vtkDataSet.h"
#include "vtkFieldData.h"
#include "vtkInformation.h"
#include "vtkInformationVector.h"
#include "vtkObjectFactory.h"
#include "vtkPolyData.h"
#include "vtkType.h"
#include "vtkUnstructuredGrid.h"
vtkStandardNewMacro(vtkDataObjectTreeToPointSetFilter);
//----------------------------------------------------------------------------
vtkDataObjectTreeToPointSetFilter::vtkDataObjectTreeToPointSetFilter()
: SubTreeCompositeIndex(0)
, MergePoints(true)
, Tolerance(0.0)
, ToleranceIsAbsolute(false)
, OutputDataSetType(VTK_UNSTRUCTURED_GRID)
{
}
//----------------------------------------------------------------------------
vtkDataObjectTreeToPointSetFilter::~vtkDataObjectTreeToPointSetFilter()
{
}
//----------------------------------------------------------------------------
int vtkDataObjectTreeToPointSetFilter::RequestDataObject(vtkInformation* vtkNotUsed(request),
vtkInformationVector** inputVector, vtkInformationVector* outputVector)
{
vtkInformation* inInfo = inputVector[0]->GetInformationObject(0);
if (!inInfo)
{
return 0;
}
if (this->OutputDataSetType != VTK_POLY_DATA && this->OutputDataSetType != VTK_UNSTRUCTURED_GRID)
{
vtkErrorMacro(
"Output type '" << vtkDataObjectTypes::GetClassNameFromTypeId(this->OutputDataSetType)
<< "' is not supported.");
return 0;
}
auto input = vtkDataObjectTree::SafeDownCast(inInfo->Get(vtkDataObject::DATA_OBJECT()));
if (input)
{
vtkInformation* info = outputVector->GetInformationObject(0);
vtkDataObject* output = info->Get(vtkDataObject::DATA_OBJECT());
if (!output || (vtkDataObjectTypes::GetTypeIdFromClassName(output->GetClassName()) !=
this->OutputDataSetType))
{
vtkSmartPointer<vtkDataObject> newOutput;
newOutput.TakeReference(vtkDataObjectTypes::NewDataObject(this->OutputDataSetType));
info->Set(vtkDataObject::DATA_OBJECT(), newOutput);
this->GetOutputPortInformation(0)->Set(
vtkDataObject::DATA_EXTENT_TYPE(), newOutput->GetExtentType());
}
return 1;
}
return 0;
}
//----------------------------------------------------------------------------
int vtkDataObjectTreeToPointSetFilter::RequestData(vtkInformation* vtkNotUsed(request),
vtkInformationVector** inputVector, vtkInformationVector* outputVector)
{
// input
vtkDataSet* ds = vtkDataSet::GetData(inputVector[0], 0);
vtkDataObjectTree* dot = vtkDataObjectTree::GetData(inputVector[0], 0);
// output
vtkDataSet* output = vtkDataSet::GetData(outputVector, 0);
vtkNew<vtkAppendDataSets> appender;
appender->SetOutputDataSetType(this->OutputDataSetType);
appender->SetMergePoints(this->MergePoints ? 1 : 0);
if (this->MergePoints)
{
appender->SetTolerance(this->Tolerance);
appender->SetToleranceIsAbsolute(this->ToleranceIsAbsolute);
}
if (ds)
{
appender->AddInputData(ds);
}
else if (dot)
{
if (this->SubTreeCompositeIndex == 0)
{
this->ExecuteSubTree(dot, appender);
}
vtkDataObjectTreeIterator* iter =
vtkDataObjectTreeIterator::SafeDownCast(dot->NewTreeIterator());
iter->VisitOnlyLeavesOff();
for (iter->InitTraversal();
!iter->IsDoneWithTraversal() && iter->GetCurrentFlatIndex() <= this->SubTreeCompositeIndex;
iter->GoToNextItem())
{
if (iter->GetCurrentFlatIndex() == this->SubTreeCompositeIndex)
{
vtkDataObject* curDO = iter->GetCurrentDataObject();
vtkDataObjectTree* curDOT = vtkDataObjectTree::SafeDownCast(curDO);
vtkDataSet* curDS = vtkUnstructuredGrid::SafeDownCast(curDO);
if (curDS && curDS->GetNumberOfPoints() > 0)
{
appender->AddInputData(curDS);
}
else if (curDOT)
{
this->ExecuteSubTree(curDOT, appender);
}
break;
}
}
iter->Delete();
}
if (appender->GetNumberOfInputConnections(0) > 0)
{
appender->Update();
output->ShallowCopy(appender->GetOutput());
// this will override field data the vtkAppendFilter passed from the first
// block. It seems like a reasonable approach, if global field data is
// present.
output->GetFieldData()->PassData(dot->GetFieldData());
}
this->RemovePartialArrays(output);
return 1;
}
//----------------------------------------------------------------------------
int vtkDataObjectTreeToPointSetFilter::FillInputPortInformation(
int vtkNotUsed(port), vtkInformation* info)
{
info->Remove(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE());
info->Append(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkDataObjectTree");
info->Append(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkDataSet");
return 1;
}
//----------------------------------------------------------------------------
void vtkDataObjectTreeToPointSetFilter::RemovePartialArrays(vtkDataSet* data)
{
vtkNew<vtkCleanArrays> cleaner;
cleaner->SetInputData(data);
cleaner->Update();
data->ShallowCopy(cleaner->GetOutput());
}
//----------------------------------------------------------------------------
void vtkDataObjectTreeToPointSetFilter::ExecuteSubTree(
vtkDataObjectTree* dot, vtkAppendDataSets* appender)
{
vtkDataObjectTreeIterator* iter = dot->NewTreeIterator();
for (iter->InitTraversal(); !iter->IsDoneWithTraversal(); iter->GoToNextItem())
{
vtkDataSet* curDS = vtkDataSet::SafeDownCast(iter->GetCurrentDataObject());
if (curDS)
{
appender->AddInputData(curDS);
}
}
iter->Delete();
}
//----------------------------------------------------------------------------
void vtkDataObjectTreeToPointSetFilter::PrintSelf(ostream& os, vtkIndent indent)
{
this->Superclass::PrintSelf(os, indent);
os << indent << "SubTreeCompositeIndex: " << this->SubTreeCompositeIndex << "\n";
os << indent << "MergePoints: " << (this->MergePoints ? "On" : "Off") << "\n";
os << indent << "Tolerance: " << this->Tolerance << "\n";
os << indent << "ToleranceIsAbsolute: " << this->ToleranceIsAbsolute << "\n";
os << indent
<< "OutputDataSetType: " << vtkDataObjectTypes::GetClassNameFromTypeId(this->OutputDataSetType)
<< "\n";
}
/*=========================================================================
Program: ParaView
Module: vtkDataObjectTreeToPointSetFilter.h
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.paraview.org/HTML/Copyright.html 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.
=========================================================================*/
/**
* @class vtkDataObjectTreeToPointSetFilter
* @brief Filter that merges blocks from a data object tree to a vtkPolyData of vtkUnstructuredGrid
*
* This filter takes a vtkDataObjectTree as input and merges the blocks in the tree into a single
* output vtkUnstructuredGrid if the ExtractSurfaces option is off. If the ExtractSurfaces option
* is on, surfaces of non-vtkPolyData datasets in the input tree will be extracted and merged to
* the output with vtkPolyData datasets.
*/
#ifndef vtkDataObjectTreeToPointSetFilter_h
#define vtkDataObjectTreeToPointSetFilter_h
#include "vtkPVVTKExtensionsRenderingModule.h" // For export macro
#include "vtkPointSetAlgorithm.h"
class vtkAppendDataSets;
class vtkDataObjectTree;
class vtkDataSet;
class VTKPVVTKEXTENSIONSRENDERING_EXPORT vtkDataObjectTreeToPointSetFilter
: public vtkPointSetAlgorithm
{
public:
static vtkDataObjectTreeToPointSetFilter* New();
vtkTypeMacro(vtkDataObjectTreeToPointSetFilter, vtkPointSetAlgorithm);
void PrintSelf(ostream& os, vtkIndent indent) override;
//@{
/**
* Get/Set the composite index of the subtree to be merged. By default set to
* 0 i.e. root, hence entire input composite dataset is merged.
*/
vtkSetMacro(SubTreeCompositeIndex, unsigned int);
vtkGetMacro(SubTreeCompositeIndex, unsigned int);
//@}
//@{
/**
* Turn on/off merging of coincidental points. Frontend to
* vtkAppendFilter::MergePoints. Default is on.
*/
vtkSetMacro(MergePoints, bool);
vtkGetMacro(MergePoints, bool);
vtkBooleanMacro(MergePoints, bool);
//@}
//@{
/**
* Get/Set the tolerance to use to find coincident points when `MergePoints`
* is `true`. Default is 0.0.
*
* This is simply passed on to the internal vtkAppendFilter::vtkLocator used to merge points.
* @sa `vtkLocator::SetTolerance`.
*/
vtkSetClampMacro(Tolerance, double, 0.0, VTK_DOUBLE_MAX);
vtkGetMacro(Tolerance, double);
//@}
//@{
/**
* Get/Set whether Tolerance is treated as an absolute or relative tolerance.
* The default is to treat it as an absolute tolerance.
*/
vtkSetMacro(ToleranceIsAbsolute, bool);
vtkGetMacro(ToleranceIsAbsolute, bool);
vtkBooleanMacro(ToleranceIsAbsolute, bool);
//@}
//@{
/**
* Get/Set the output type produced by this filter. Only blocks compatible with the output type
* will be merged in the output. For example, if the output type is vtkPolyData, then
* blocks of type vtkImageData, vtkStructuredGrid, etc. will not be merged - only vtkPolyData
* can be merged into a vtkPolyData. On the other hand, if the output type is
* vtkUnstructuredGrid, then blocks of almost any type will be merged in the output.
* Valid values are VTK_POLY_DATA and VTK_UNSTRUCTURED_GRID defined in vtkType.h.
* Defaults to VTK_UNSTRUCTURED_GRID.
*/
//@}
vtkSetMacro(OutputDataSetType, int);
vtkGetMacro(OutputDataSetType, int);
protected:
vtkDataObjectTreeToPointSetFilter();
~vtkDataObjectTreeToPointSetFilter() override;
int RequestDataObject(vtkInformation* request, vtkInformationVector** inputVector,
vtkInformationVector* outputVector) override;
int RequestData(vtkInformation* request, vtkInformationVector** inputVector,
vtkInformationVector* outputVector) override;
int FillInputPortInformation(int port, vtkInformation* info) override;
/**
* Remove point/cell arrays not present on all processes.
*/
void RemovePartialArrays(vtkDataSet* data);
void ExecuteSubTree(vtkDataObjectTree* dot, vtkAppendDataSets* appender);
unsigned int SubTreeCompositeIndex;
bool MergePoints;
double Tolerance;
bool ToleranceIsAbsolute;
int OutputDataSetType;
private:
vtkDataObjectTreeToPointSetFilter(const vtkDataObjectTreeToPointSetFilter&) = delete;
void operator=(const vtkDataObjectTreeToPointSetFilter&) = delete;
};
#endif
Markdown is supported
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