Commit db3c55e9 authored by Will Schroeder's avatar Will Schroeder Committed by Kitware Robot

Merge topic 'GeneralizedPointInterpolator'

4173dbea Proper pipeline behavior
045eac98 Test cleanup
956b5591 Refactored kernel API to support different basis
eadf3ec1 Created new Filters/Points module
ae3d4315 Preparation for reorganization
e5c4dd6d Added new class for oriented gaussian distributions
b92e2cbc Compiler warnings
0fb06f89 Compiler warnings
...
Acked-by: Kitware Robot's avatarKitware Robot <kwrobot@kitware.com>
Reviewed-by: Berk Geveci's avatarBerk Geveci <berk.geveci@kitware.com>
Merge-request: !1062
parents 9e4cfa2b 4173dbea
......@@ -36,6 +36,7 @@ vtk_add_test_cxx(${vtk-module}CxxTests tests
TestTreeBFSIterator.cxx
TestTreeDFSIterator.cxx
TestTriangle.cxx
TimePointLocators.cxx
otherCellArray.cxx
otherCellBoundaries.cxx
otherCellPosition.cxx
......
/*=========================================================================
Program: Visualization Toolkit
Module: TestPointLocators.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 "vtkPointLocator.h"
#include "vtkStaticPointLocator.h"
#include "vtkKdTree.h"
#include "vtkKdTreePointLocator.h"
#include "vtkOctreePointLocator.h"
#include "vtkPoints.h"
#include "vtkPolyData.h"
#include "vtkTimerLog.h"
#include "vtkMath.h"
int TimePointLocators(int , char *[])
{
int nPts = 1000000;
int nQ = nPts/10;
int N = 10;
double R = 0.01;
vtkTimerLog* timer = vtkTimerLog::New();
double buildTime[4], cpTime[4], cnpTime[4], crpTime[4];
for (int i=0; i<4; ++i)
{
buildTime[i] = cpTime[i] = cnpTime[i] = crpTime[i] = 0.0;
}
cout << "\nTiming for " << nPts << " points, " << nQ << " queries\n";
// Populate a list of points and query locations
vtkPoints* points = vtkPoints::New();
points->SetDataTypeToDouble();
points->SetNumberOfPoints(nPts);
for (int i=0; i<nPts; ++i)
{
points->SetPoint(i, vtkMath::Random(-1,1), vtkMath::Random(-1,1),
vtkMath::Random(-1,1));
}
vtkPolyData* polydata = vtkPolyData::New();
polydata->SetPoints(points);
points->ComputeBounds();
// Create points array which are positions to probe data with FindClosestPoint().
vtkPoints* qPoints = vtkPoints::New();
qPoints->SetDataTypeToDouble();
qPoints->SetNumberOfPoints(nQ);
vtkMath::RandomSeed(314159);
for (int i=0; i<nQ; ++i)
{
qPoints->SetPoint(i, vtkMath::Random(-1,1), vtkMath::Random(-1,1),
vtkMath::Random(-1,1));
}
vtkIdList *closest = vtkIdList::New();
//---------------------------------------------------------------------------
// The simple uniform binning point locator
vtkPointLocator* uniformLocator = vtkPointLocator::New();
timer->StartTimer();
uniformLocator->SetDataSet(polydata);
uniformLocator->BuildLocator();
timer->StopTimer();
uniformLocator->Delete();
buildTime[0] = timer->GetElapsedTime();
uniformLocator = vtkPointLocator::New();
uniformLocator->SetDataSet(polydata);
uniformLocator->BuildLocator();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
uniformLocator->FindClosestPoint(qPoints->GetPoint(i));
}
timer->StopTimer();
cpTime[0] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
uniformLocator->FindClosestNPoints(N, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
cnpTime[0] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
uniformLocator->FindPointsWithinRadius(R, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
crpTime[0] = timer->GetElapsedTime();
uniformLocator->Delete();
//---------------------------------------------------------------------------
// The simple uniform binning point locator, this time built
// statically and threaded (may be much faster on threaded machine).
vtkStaticPointLocator* staticLocator = vtkStaticPointLocator::New();
timer->StartTimer();
staticLocator->SetDataSet(polydata);
staticLocator->BuildLocator();
timer->StopTimer();
staticLocator->Delete();
buildTime[1] = timer->GetElapsedTime();
staticLocator = vtkStaticPointLocator::New();
staticLocator->SetDataSet(polydata);
staticLocator->BuildLocator();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
staticLocator->FindClosestPoint(qPoints->GetPoint(i));
}
timer->StopTimer();
cpTime[1] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
staticLocator->FindClosestNPoints(N, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
cnpTime[1] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
staticLocator->FindPointsWithinRadius(R, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
crpTime[1] = timer->GetElapsedTime();
staticLocator->Delete();
//---------------------------------------------------------------------------
// The KD Tree point locator
vtkKdTreePointLocator* kdTreeLocator = vtkKdTreePointLocator::New();
timer->StartTimer();
kdTreeLocator->SetDataSet(polydata);
kdTreeLocator->BuildLocator();
kdTreeLocator->Delete();
timer->StopTimer();
buildTime[2] = timer->GetElapsedTime();
kdTreeLocator = vtkKdTreePointLocator::New();
kdTreeLocator->SetDataSet(polydata);
kdTreeLocator->BuildLocator();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
kdTreeLocator->FindClosestPoint(qPoints->GetPoint(i));
}
timer->StopTimer();
cpTime[2] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
kdTreeLocator->FindClosestNPoints(N, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
cnpTime[2] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
kdTreeLocator->FindPointsWithinRadius(R, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
crpTime[2] = timer->GetElapsedTime();
kdTreeLocator->Delete();
//---------------------------------------------------------------------------
// Octree point locator
vtkOctreePointLocator* octreeLocator = vtkOctreePointLocator::New();
timer->StartTimer();
octreeLocator->SetDataSet(polydata);
octreeLocator->BuildLocator();
octreeLocator->Delete();
timer->StopTimer();
buildTime[3] = timer->GetElapsedTime();
octreeLocator = vtkOctreePointLocator::New();
octreeLocator->SetDataSet(polydata);
octreeLocator->BuildLocator();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
octreeLocator->FindClosestPoint(qPoints->GetPoint(i));
}
timer->StopTimer();
cpTime[3] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
octreeLocator->FindClosestNPoints(N, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
cnpTime[3] = timer->GetElapsedTime();
timer->StartTimer();
for (int i=0; i<nQ; ++i)
{
octreeLocator->FindPointsWithinRadius(R, qPoints->GetPoint(i), closest);
}
timer->StopTimer();
crpTime[3] = timer->GetElapsedTime();
octreeLocator->Delete();
//---------------------------------------------------------------------------
// Print out the statistics
cout << "Build and delete tree\n";
cout << "\tUniform: " << buildTime[0] << "\n";
cout << "\tStatic: " << buildTime[1] << "\n";
cout << "\tOctree: " << buildTime[2] << "\n";
cout << "\tKD Tree: " << buildTime[3] << "\n";
cout << "Closest point queries\n";
cout << "\tUniform: " << cpTime[0] << "\n";
cout << "\tStatic: " << cpTime[1] << "\n";
cout << "\tOctree: " << cpTime[2] << "\n";
cout << "\tKD Tree: " << cpTime[3] << "\n";
cout << "Closest N points queries\n";
cout << "\tUniform: " << cnpTime[0] << "\n";
cout << "\tStatic: " << cnpTime[1] << "\n";
cout << "\tOctree: " << cnpTime[2] << "\n";
cout << "\tKD Tree: " << cnpTime[3] << "\n";
cout << "Closest points within radius queries\n";
cout << "\tUniform: " << crpTime[0] << "\n";
cout << "\tStatic: " << crpTime[1] << "\n";
cout << "\tOctree: " << crpTime[2] << "\n";
cout << "\tKD Tree: " << crpTime[3] << "\n";
cout << "Total time\n";
cout << "\tUniform: " << (buildTime[0] + cpTime[0] + cnpTime[0] + crpTime[0]) << "\n";
cout << "\tStatic: " << (buildTime[1] + cpTime[1] + cnpTime[1] + crpTime[1]) << "\n";
cout << "\tOctree: " << (buildTime[2] + cpTime[2] + cnpTime[2] + crpTime[2]) << "\n";
cout << "\tKD Tree: " << (buildTime[3] + cpTime[3] + cnpTime[3] + crpTime[3]) << "\n";
timer->Delete();
points->Delete();
qPoints->Delete();
polydata->Delete();
closest->Delete();
// Always return success, although the test infrastructure should catch
// excessive execution times.
return 0;
}
/*=========================================================================
Program: Visualization Toolkit
Module: vtkVoronoiKernel.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 "vtkVoronoiKernel.h"
#include "vtkAbstractPointLocator.h"
//----------------------------------------------------------------------------
vtkVoronoiKernel::vtkVoronoiKernel()
{
}
//----------------------------------------------------------------------------
vtkVoronoiKernel::~vtkVoronoiKernel()
{
}
//----------------------------------------------------------------------------
void vtkVoronoiKernel::
vtkIdType ComputeWeights(double x[3], vtkIdList *pIds,vtkDoubleArray *weights)
{
}
//----------------------------------------------------------------------------
void vtkVoronoiKernel::PrintSelf(ostream& os, vtkIndent indent)
{
this->Superclass::PrintSelf(os,indent);
}
set(Module_SRCS
vtkEllipsoidalGaussianKernel.cxx
vtkGaussianKernel.cxx
vtkInterpolationKernel.cxx
vtkPointInterpolator.cxx
vtkLinearKernel.cxx
vtkSPHKernel.cxx
vtkShepardKernel.cxx
vtkVoronoiKernel.cxx
)
set_source_files_properties(
vtkInterpolationKernel
ABSTRACT
)
vtk_module_library(vtkFiltersPoints ${Module_SRCS})
vtk_add_test_python(
TestPointInterpolator.py
TestPointInterpolator2.py
)
#!/usr/bin/env python
import vtk
from vtk.test import Testing
from vtk.util.misc import vtkGetDataRoot
VTK_DATA_ROOT = vtkGetDataRoot()
# Parameters for debugging
res = 1000
# create pipeline
#
pl3d = vtk.vtkMultiBlockPLOT3DReader()
pl3d.SetXYZFileName(VTK_DATA_ROOT + "/Data/combxyz.bin")
pl3d.SetQFileName(VTK_DATA_ROOT + "/Data/combq.bin")
pl3d.SetScalarFunctionNumber(100)
pl3d.SetVectorFunctionNumber(202)
pl3d.Update()
output = pl3d.GetOutput().GetBlock(0)
# Create a probe plane
center = output.GetCenter()
plane = vtk.vtkPlaneSource()
plane.SetResolution(res,res)
plane.SetOrigin(0,0,0)
plane.SetPoint1(10,0,0)
plane.SetPoint2(0,10,0)
plane.SetCenter(center)
plane.SetNormal(0,1,0)
# Reuse the locator
locator = vtk.vtkStaticPointLocator()
locator.SetDataSet(output)
locator.BuildLocator()
# Voronoi kernel------------------------------------------------
voronoiKernel = vtk.vtkVoronoiKernel()
interpolator = vtk.vtkPointInterpolator()
interpolator.SetInputConnection(plane.GetOutputPort())
interpolator.SetSourceData(output)
interpolator.SetKernel(voronoiKernel)
interpolator.SetLocator(locator)
# Time execution
timer = vtk.vtkTimerLog()
timer.StartTimer()
interpolator.Update()
timer.StopTimer()
time = timer.GetElapsedTime()
print("Interpolate Points (Voronoi): {0}".format(time))
intMapper = vtk.vtkPolyDataMapper()
intMapper.SetInputConnection(interpolator.GetOutputPort())
intActor = vtk.vtkActor()
intActor.SetMapper(intMapper)
# Create an outline
outline = vtk.vtkStructuredGridOutlineFilter()
outline.SetInputData(output)
outlineMapper = vtk.vtkPolyDataMapper()
outlineMapper.SetInputConnection(outline.GetOutputPort())
outlineActor = vtk.vtkActor()
outlineActor.SetMapper(outlineMapper)
# Gaussian kernel-------------------------------------------------------
gaussianKernel = vtk.vtkGaussianKernel()
#gaussianKernel = vtk.vtkEllipsoidalGaussianKernel()
#gaussianKernel.UseScalarsOn()
#gaussianKernel.UseNormalsOn()
gaussianKernel.SetSharpness(4)
gaussianKernel.SetRadius(0.5)
interpolator1 = vtk.vtkPointInterpolator()
interpolator1.SetInputConnection(plane.GetOutputPort())
interpolator1.SetSourceData(output)
interpolator1.SetKernel(gaussianKernel)
interpolator1.SetLocator(locator)
interpolator1.SetNullPointsStrategyToNullValue()
# Time execution
timer.StartTimer()
interpolator1.Update()
timer.StopTimer()
time = timer.GetElapsedTime()
print("Interpolate Points (Gaussian): {0}".format(time))
intMapper1 = vtk.vtkPolyDataMapper()
intMapper1.SetInputConnection(interpolator1.GetOutputPort())
intActor1 = vtk.vtkActor()
intActor1.SetMapper(intMapper1)
# Create an outline
outline1 = vtk.vtkStructuredGridOutlineFilter()
outline1.SetInputData(output)
outlineMapper1 = vtk.vtkPolyDataMapper()
outlineMapper1.SetInputConnection(outline1.GetOutputPort())
outlineActor1 = vtk.vtkActor()
outlineActor1.SetMapper(outlineMapper1)
# Shepard kernel-------------------------------------------------------
shepardKernel = vtk.vtkShepardKernel()
shepardKernel.SetPowerParameter(2)
shepardKernel.SetRadius(0.5)
interpolator2 = vtk.vtkPointInterpolator()
interpolator2.SetInputConnection(plane.GetOutputPort())
interpolator2.SetSourceData(output)
interpolator2.SetKernel(shepardKernel)
interpolator2.SetLocator(locator)
interpolator2.SetNullPointsStrategyToMaskPoints()
# Time execution
timer.StartTimer()
interpolator2.Update()
timer.StopTimer()
time = timer.GetElapsedTime()
print("Interpolate Points (Shepard): {0}".format(time))
intMapper2 = vtk.vtkPolyDataMapper()
intMapper2.SetInputConnection(interpolator2.GetOutputPort())
intActor2 = vtk.vtkActor()
intActor2.SetMapper(intMapper2)
# Create an outline
outline2 = vtk.vtkStructuredGridOutlineFilter()
outline2.SetInputData(output)
outlineMapper2 = vtk.vtkPolyDataMapper()
outlineMapper2.SetInputConnection(outline2.GetOutputPort())
outlineActor2 = vtk.vtkActor()
outlineActor2.SetMapper(outlineMapper2)
# SPH kernel-------------------------------------------------------
SPHKernel = vtk.vtkSPHKernel()
interpolator3 = vtk.vtkPointInterpolator()
interpolator3.SetInputConnection(plane.GetOutputPort())
interpolator3.SetSourceData(output)
interpolator3.SetKernel(SPHKernel)
interpolator3.SetLocator(locator)
# Time execution
timer.StartTimer()
interpolator3.Update()
timer.StopTimer()
time = timer.GetElapsedTime()
print("Interpolate Points (SPH): {0}".format(time))
intMapper3 = vtk.vtkPolyDataMapper()
intMapper3.SetInputConnection(interpolator3.GetOutputPort())
intActor3 = vtk.vtkActor()
intActor3.SetMapper(intMapper3)
# Create an outline
outline3 = vtk.vtkStructuredGridOutlineFilter()
outline3.SetInputData(output)
outlineMapper3 = vtk.vtkPolyDataMapper()
outlineMapper3.SetInputConnection(outline3.GetOutputPort())
outlineActor3 = vtk.vtkActor()
outlineActor3.SetMapper(outlineMapper3)
# Create the RenderWindow, Renderer and both Actors
#
ren0 = vtk.vtkRenderer()
ren0.SetViewport(0,0,.5,.5)
ren1 = vtk.vtkRenderer()
ren1.SetViewport(0.5,0,1,.5)
ren2 = vtk.vtkRenderer()
ren2.SetViewport(0,0.5,.5,1)
ren3 = vtk.vtkRenderer()
ren3.SetViewport(0.5,0.5,1,1)
renWin = vtk.vtkRenderWindow()
renWin.AddRenderer(ren0)
renWin.AddRenderer(ren1)
renWin.AddRenderer(ren2)
renWin.AddRenderer(ren3)
iren = vtk.vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)
# Add the actors to the renderer, set the background and size
#
ren0.AddActor(intActor)
ren0.AddActor(outlineActor)
ren0.SetBackground(0.1, 0.2, 0.4)
ren1.AddActor(intActor1)
ren1.AddActor(outlineActor1)
ren1.SetBackground(0.1, 0.2, 0.4)
ren2.AddActor(intActor2)
ren2.AddActor(outlineActor2)
ren2.SetBackground(0.1, 0.2, 0.4)
ren3.AddActor(intActor3)
ren3.AddActor(outlineActor3)
ren3.SetBackground(0.1, 0.2, 0.4)
renWin.SetSize(500, 500)
cam = ren0.GetActiveCamera()
cam.SetClippingRange(3.95297, 50)
cam.SetFocalPoint(8.88908, 0.595038, 29.3342)
cam.SetPosition(-12.3332, 31.7479, 41.2387)
cam.SetViewUp(0.060772, -0.319905, 0.945498)
ren1.SetActiveCamera(cam)
ren2.SetActiveCamera(cam)
ren3.SetActiveCamera(cam)
iren.Initialize()
# render the image
#
renWin.Render()
#iren.Start()
#!/usr/bin/env python
import vtk
from vtk.test import Testing
from vtk.util.misc import vtkGetDataRoot
VTK_DATA_ROOT = vtkGetDataRoot()
# Interpolate onto a volume
# Parameters for debugging
res = 100
# create pipeline
#
extent = [0,56, 0,32, 0,24]
pl3d = vtk.vtkMultiBlockPLOT3DReader()
pl3d.SetXYZFileName(VTK_DATA_ROOT + "/Data/combxyz.bin")
pl3d.SetQFileName(VTK_DATA_ROOT + "/Data/combq.bin")
pl3d.SetScalarFunctionNumber(100)
pl3d.SetVectorFunctionNumber(202)
pl3d.Update()
output = pl3d.GetOutput().GetBlock(0)
# Create a probe volume
center = output.GetCenter()
bounds = output.GetBounds()
length = output.GetLength()
probe = vtk.vtkImageData()
probe.SetDimensions(res,res,res)
probe.SetOrigin(bounds[0],bounds[2],bounds[4])
probe.SetSpacing((bounds[1]-bounds[0])/(res-1),
(bounds[3]-bounds[2])/(res-1),
(bounds[5]-bounds[4])/(res-1))
# Reuse the locator
locator = vtk.vtkStaticPointLocator()
locator.SetDataSet(output)
locator.BuildLocator()
# Use a gaussian kernel------------------------------------------------
gaussianKernel = vtk.vtkGaussianKernel()
gaussianKernel.SetRadius(0.5)
gaussianKernel.SetSharpness(4)
print ("Radius: {0}".format(gaussianKernel.GetRadius()))
interpolator = vtk.vtkPointInterpolator()
interpolator.SetInputData(probe)
interpolator.SetSourceData(output)
interpolator.SetKernel(gaussianKernel)
interpolator.SetLocator(locator)
# Time execution
timer = vtk.vtkTimerLog()
timer.StartTimer()
interpolator.Update()
timer.StopTimer()
time = timer.GetElapsedTime()
print("Interpolate Points (Volume probe): {0}".format(time))
intMapper = vtk.vtkDataSetMapper()
intMapper.SetInputConnection(interpolator.GetOutputPort())
intActor = vtk.vtkActor()
intActor.SetMapper(intMapper)
# Create an outline
outline = vtk.vtkStructuredGridOutlineFilter()
outline.SetInputData(output)
outlineMapper = vtk.vtkPolyDataMapper()
outlineMapper.SetInputConnection(outline.GetOutputPort())
outlineActor = vtk.vtkActor()
outlineActor.SetMapper(outlineMapper)
# Create the RenderWindow, Renderer and both Actors
#
ren0 = vtk.vtkRenderer()
renWin = vtk.vtkRenderWindow()
renWin.AddRenderer(ren0)
iren = vtk.vtkRenderWindowInteractor()
iren.SetRenderWindow(renWin)
# Add the actors to the renderer, set the background and size
#
ren0.AddActor(intActor)
ren0.AddActor(outlineActor)
ren0.SetBackground(0.1, 0.2, 0.4)
renWin.SetSize(250,250)
cam = ren0.GetActiveCamera()
cam.SetClippingRange(3.95297, 50)
cam.SetFocalPoint(8.88908, 0.595038, 29.3342)
cam.SetPosition(-12.3332, 31.7479, 41.2387)
cam.SetViewUp(0.060772, -0.319905, 0.945498)
iren.Initialize()
# render the image
#
renWin.Render()
#iren.Start()
vtk_module(vtkFiltersPoints
GROUPS
StandAlone
DEPENDS
vtkCommonExecutionModel
vtkCommonSystem
vtkCommonMisc
vtkCommonTransforms
vtkCommonMath