Commit 03de2912 authored by T.J. Corona's avatar T.J. Corona Committed by Kitware Robot
Browse files

Merge topic 'unstructured-grid-quadric-decimation'

4665a1b9

 Add vtkUnstructuredGridQuadricDecimation.
Acked-by: Kitware Robot's avatarKitware Robot <kwrobot@kitware.com>
Merge-request: !1600
parents 4e546820 4665a1b9
Pipeline #20725 failed with stage
in 0 seconds
......@@ -64,6 +64,7 @@ set(Module_SRCS
vtkTransposeTable.cxx
vtkTriangleFilter.cxx
vtkTubeFilter.cxx
vtkUnstructuredGridQuadricDecimation.cxx
vtkVectorDot.cxx
vtkVectorNorm.cxx
vtkWindowedSincPolyDataFilter.cxx
......
......@@ -38,6 +38,7 @@ vtk_add_test_cxx(${vtk-module}CxxTests tests
TestThresholdPoints.cxx,NO_VALID
TestTransposeTable.cxx,NO_VALID
TestTubeFilter.cxx,NO_VALID
TestUnstructuredGridQuadricDecimation.cxx,NO_VALID
UnitTestMaskPoints.cxx,NO_VALID
UnitTestMergeFilter.cxx,NO_VALID
)
......
/*=========================================================================
Program: Visualization Toolkit
Module: TestUnstructuredGridQuadricDecimation.cxx
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.kitware.com/Copyright.htm 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 <vtkActor.h>
#include <vtkCellArray.h>
#include <vtkCleanPolyData.h>
#include <vtkDataSetMapper.h>
#include <vtkDelaunay3D.h>
#include <vtkDoubleArray.h>
#include <vtkMath.h>
#include <vtkPointData.h>
#include <vtkPoints.h>
#include <vtkPointSource.h>
#include <vtkPolyData.h>
#include <vtkPolygon.h>
#include <vtkProperty.h>
#include <vtkRenderer.h>
#include <vtkRenderWindow.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkSmartPointer.h>
#include <vtkUnstructuredGrid.h>
#include <vtkUnstructuredGridQuadricDecimation.h>
int TestUnstructuredGridQuadricDecimation(int, char *[])
{
// This test constructs a tetrahedrally meshed sphere by first generating
// <numberOfOriginalPoints> points randomly placed within a unit sphere, then
// removing points that overlap within a tolerance, and finally constructing a
// delaunay 3d tetrahedralization from the points. Additionally, point data
// corresponding to the points distance from the origin are added to this
// data. The resulting tetrahedral mesh is then decimated <numberOfTests>
// times, each time with a target reduction facter <targetReduction[test]>.
// The number of remaining tetrahedra is then compared to the original number
// of tetrahedra and compared against the target reduction factor. If the
// difference is greater than <absTolerance>, the test fails. Otherwise, the
// test passes.
// # of points to generate the original tetrahedral mesh
const vtkIdType numberOfOriginalPoints = 1.e4;
// # of decimation tests to perform
const vtkIdType numberOfTests = 4;
// target reduction values for each test
const double targetReduction[numberOfTests] = {.1,.3,.5,.7};
// absolute tolerance between the expected and received tetrahedron reduction
// to determine whether the decimation successfully executed
const double absTolerance = 1.e-1;
// Generate points within a unit sphere centered at the origin.
vtkSmartPointer<vtkPointSource> source =
vtkSmartPointer<vtkPointSource>::New();
source->SetNumberOfPoints(numberOfOriginalPoints);
source->SetCenter(0.,0.,0.);
source->SetRadius(1.);
source->SetDistributionToUniform();
source->SetOutputPointsPrecision(vtkAlgorithm::DOUBLE_PRECISION);
// Clean the polydata. This will remove overlapping points that may be
// present in the input data.
vtkSmartPointer<vtkCleanPolyData> cleaner =
vtkSmartPointer<vtkCleanPolyData>::New();
cleaner->SetInputConnection(source->GetOutputPort());
cleaner->Update();
// Create point data for use in decimation (the point data acts as a fourth
// dimension in a Euclidean metric for determining the "nearness" of points).
vtkPolyData* pd = cleaner->GetOutput();
vtkPoints* points = pd->GetPoints();
vtkSmartPointer<vtkDoubleArray> radius =
vtkSmartPointer<vtkDoubleArray>::New();
radius->SetName("radius");
radius->SetNumberOfComponents(1);
radius->SetNumberOfTuples(points->GetNumberOfPoints());
double xyz[3];
double r;
for (vtkIdType i=0;i<points->GetNumberOfPoints();i++)
{
points->GetPoint(i,xyz);
r = std::sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1] + xyz[2]*xyz[2]);
radius->SetTypedTuple(i,&r);
}
pd->GetPointData()->SetScalars(radius);
// Generate a tetrahedral mesh from the input points. By
// default, the generated volume is the convex hull of the points.
vtkSmartPointer<vtkDelaunay3D> delaunay3D =
vtkSmartPointer<vtkDelaunay3D>::New();
delaunay3D->SetInputData(pd);
delaunay3D->Update();
const vtkIdType numberOfOriginalTetras =
delaunay3D->GetOutput()->GetNumberOfCells();
for (vtkIdType test = 0; test < numberOfTests; test++)
{
// Decimate the tetrahedral mesh.
vtkSmartPointer<vtkUnstructuredGridQuadricDecimation> decimate =
vtkSmartPointer<vtkUnstructuredGridQuadricDecimation>::New();
decimate->SetInputConnection(delaunay3D->GetOutputPort());
decimate->SetScalarsName("radius");
decimate->SetTargetReduction(targetReduction[test]);
decimate->Update();
// Compare the resultant decimation fraction with the expected fraction.
double fraction =
(1. - static_cast<double>(decimate->GetOutput()
->GetNumberOfCells()) / numberOfOriginalTetras);
std::cout<<"Test # "<<test<<std::endl;
std::cout<<"number of original tetras: "<<numberOfOriginalTetras<<std::endl;
std::cout<<"number of tetras after decimation: "<<decimate->GetOutput()
->GetNumberOfCells()<<std::endl;
std::cout<<"fraction: "<<fraction<<std::endl;
std::cout<<"expected fraction: "<<targetReduction[test]<<std::endl;
if (std::fabs(fraction - targetReduction[test]) > absTolerance)
{
return EXIT_FAILURE;
}
}
return EXIT_SUCCESS;
}
This diff is collapsed.
/*=========================================================================
Program: Visualization Toolkit
Module: vtkVectorNorm.h
Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
All rights reserved.
See Copyright.txt or http://www.kitware.com/Copyright.htm 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 2007, 2008 by University of Utah.
=========================================================================*/
// .NAME vtkUnstructuredGridQuadricDecimation - reduce the number of
// tetrahedra in a mesh
//
// .SECTION Description
//
// vtkUnstructuredGridQuadricDecimation is a class that simplifies
// tetrahedral meshes using randomized multiple choice edge
// collapses. The input to this filter is a vtkUnstructuredGrid object
// with a single scalar field (specifying in the ScalarsName
// attribute). Users can determine the size of the output mesh by
// either setting the value of TargetReduction or
// NumberOfTetsOutput. The BoundaryWeight can be set to control how
// well the mesh boundary should be preserved. The implementation uses
// Michael Garland's generalized Quadric Error Metrics, the Corner
// Table representation and the Standard Conjugate Gradient Method to
// order the edge collapse sequence.
//
// Instead of using the traditional priority queue, the algorithm uses
// a randomized approach to yield faster performance with comparable
// quality. At each step, a set of 8 random candidate edges are
// examined to select the best edge to collapse. This number can also
// be changed by users through NumberOfCandidates.
//
// For more information as well as the streaming version of this
// algorithm see:
//
// "Streaming Simplification of Tetrahedral Meshes" by H. T. Vo,
// S. P. Callahan, P. Lindstrom, V. Pascucci and C. T. Silva, IEEE
// Transactions on Visualization and Computer Graphics, 13(1):145-155,
// 2007.
//
// .SECTION Acknowledgments
//
// This code was developed by Huy T. Vo under the supervision of
// Prof. Claudio T. Silva. The code also contains contributions from
// Peter Lindstrom and Steven P. Callahan.
//
// The work was supported by grants, contracts, and gifts from the
// National Science Foundation, the Department of Energy and IBM.
#ifndef vtkUnstructuredGridQuadricDecimation_h
#define vtkUnstructuredGridQuadricDecimation_h
#include "vtkFiltersCoreModule.h" // For export macro
#include "vtkUnstructuredGridAlgorithm.h"
class VTKFILTERSCORE_EXPORT vtkUnstructuredGridQuadricDecimation : public vtkUnstructuredGridAlgorithm
{
public:
vtkTypeMacro(vtkUnstructuredGridQuadricDecimation, vtkUnstructuredGridAlgorithm);
void PrintSelf(ostream& os, vtkIndent indent);
static vtkUnstructuredGridQuadricDecimation *New();
// The following 3 parameters will control the process of simplification in
// the priority:
// NumberOfEdgesToDecimate, NumberOfTetsOutput, TargetReduction.
// If NumberOfEdgesToDecimate is 0, NumberOfTetsOutput will be considered. If
// NumbersOfTetsOutput is also 0, then TargetReduction will control the
// output.
// Description:
// Set/Get the desired reduction (expressed as a fraction of the original
// number of tetrehedra)
vtkSetMacro(TargetReduction, double);
vtkGetMacro(TargetReduction, double);
// Description:
// Set/Get the desired number of tetrahedra to be outputed
vtkSetMacro(NumberOfTetsOutput, int);
vtkGetMacro(NumberOfTetsOutput, int);
// Description:
// Set/Get the desired number of edge to collapse
vtkSetMacro(NumberOfEdgesToDecimate, int);
vtkGetMacro(NumberOfEdgesToDecimate, int);
// Description:
// Set/Get the number of candidates selected for each randomized set before
// performing an edge collapse. Increasing this number can help producing
// higher quality output but it will be slower. Default is 8.
vtkSetMacro(NumberOfCandidates, int);
vtkGetMacro(NumberOfCandidates, int);
// Description:
// Enable(1)/Disable(0) the feature of temporarily doubling the number of
// candidates for each randomized set if the quadric error was significantly
// increased over the last edge collapse, i.e. if the ratio between the error
// difference and the last error is over some threshold. Basically, we are
// trying to make up for cases when random selection returns so many 'bad'
// edges. By doing this we can achieve a higher quality output with much less
// time than just double the NumberOfCandidates. Default is Enabled(1)
vtkSetMacro(AutoAddCandidates, int);
vtkGetMacro(AutoAddCandidates, int);
// Description:
// Set/Get the threshold that decides when to double the set size.
// Default is 0.4.
vtkSetMacro(AutoAddCandidatesThreshold, double);
vtkGetMacro(AutoAddCandidatesThreshold, double);
// Description:
// Set/Get the weight of the boundary on the quadric metrics. The larger
// the number, the better the boundary is preserved.
vtkSetMacro(BoundaryWeight, double);
vtkGetMacro(BoundaryWeight, double);
// Description:
// Set/Get the scalar field name used for simplification
vtkSetStringMacro(ScalarsName);
vtkGetStringMacro(ScalarsName);
enum
{
NO_ERROR=0,
NON_TETRAHEDRA=1,
NO_SCALARS=2,
NO_CELLS=3
};
protected:
vtkUnstructuredGridQuadricDecimation();
~vtkUnstructuredGridQuadricDecimation();
void ReportError(int err);
int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *);
int NumberOfTetsOutput;
int NumberOfEdgesToDecimate;
int NumberOfCandidates;
int AutoAddCandidates;
double TargetReduction;
double AutoAddCandidatesThreshold;
double BoundaryWeight;
char *ScalarsName;
private:
vtkUnstructuredGridQuadricDecimation(const vtkUnstructuredGridQuadricDecimation&) VTK_DELETE_FUNCTION;
void operator=(const vtkUnstructuredGridQuadricDecimation&) VTK_DELETE_FUNCTION;
};
#endif
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