Add new points in the local point neighborhood

In the spirit of mesh subdivision filter, this filter operates on
a point cloud, generating new points in the local neighborhood of
existing points.

Common/DataModel/vtkArrayListTemplate.h

@@ -226,6 +226,10 @@ struct ArrayList
                  vtkDataSetAttributes *outPD, double nullValue=0.0,
                  bool promote=true);
 
+  // Add an array that interpolates from its own attribute values
+  void AddSelfInterpolatingArrays(vtkIdType numOutPts, vtkDataSetAttributes *attr,
+                                  double nullValue=0.0);
+
   // Add a pair of arrays (manual insertion). Returns the output array created,
   // if any. No array may be created if \c inArray was previously marked as
   // excluded using ExcludeArray().

Common/DataModel/vtkArrayListTemplate.txx

@@ -181,4 +181,37 @@ AddArrays(vtkIdType numOutPts, vtkDataSetAttributes *inPD, vtkDataSetAttributes
   }//for each candidate array
 }
 
+//----------------------------------------------------------------------------
+// Add the arrays to interpolate here. This presumes that vtkDataSetAttributes::CopyData() or
+// vtkDataSetAttributes::InterpolateData() has been called. This special version creates an
+// array pair that interpolates from itself.
+inline void ArrayList::
+AddSelfInterpolatingArrays(vtkIdType numOutPts, vtkDataSetAttributes *attr, double nullValue)
+{
+  // Build the vector of interpolation pairs. Note that CopyAllocate/InterpolateAllocate should have
+  // been called at this point (output arrays created and allocated).
+  vtkDataArray *iArray;
+  int iType, iNumComp;
+  void *iD;
+  int i, numArrays = attr->GetNumberOfArrays();
+
+  for (i=0; i < numArrays; ++i)
+  {
+    iArray = attr->GetArray(i);
+    if ( iArray && ! this->IsExcluded(iArray) )
+    {
+      iType = iArray->GetDataType();
+      iNumComp = iArray->GetNumberOfComponents();
+      iArray->WriteVoidPointer(0,numOutPts*iNumComp); //allocates memory, preserves data
+      iD = iArray->GetVoidPointer(0);
+      switch (iType)
+      {
+        vtkTemplateMacro(CreateArrayPair(this, static_cast<VTK_TT *>(iD),
+                         static_cast<VTK_TT *>(iD),numOutPts,iNumComp,
+                         iArray,static_cast<VTK_TT>(nullValue)));
+      }//over all VTK types
+    }//if not excluded
+  }//for each candidate array
+}
+
 #endif

Filters/Points/CMakeLists.txt

@@ -1,5 +1,6 @@
 set(Module_SRCS
   vtkBoundedPointSource.cxx
+  vtkDensifyPointCloudFilter.cxx
   vtkEllipsoidalGaussianKernel.cxx
   vtkEuclideanClusterExtraction.cxx
   vtkExtractHierarchicalBins.cxx

Filters/Points/Testing/Data/Baseline/TestDensifyPointCloudFilter.png.md5

@@ -0,0 +1 @@
+165f487c14fb4e86f5b75ebefc70e95e

Filters/Points/Testing/Data/Baseline/TestDensifyPointCloudFilter2.png.md5

@@ -0,0 +1 @@
+e547b0ea93a2562b66af0f5355a85fc3

Filters/Points/Testing/Python/CMakeLists.txt

@@ -1,4 +1,5 @@
 vtk_add_test_python(
+  TestDensifyPointCloudFilter2.py
   TestEllipsoidalGaussianKernel.py
   TestPointInterpolator.py
   TestPointInterpolator2.py
@@ -10,6 +11,7 @@ vtk_add_test_python(
 
 if ("${VTK_RENDERING_BACKEND}" STREQUAL "OpenGL2")
   vtk_add_test_python(
+    TestDensifyPointCloudFilter.py
     TestEuclideanClusterExtraction.py
     TestEuclideanClusterExtraction2.py
     TestExtractPoints.py

Filters/Points/Testing/Python/TestDensifyPointCloudFilter.py

@@ -0,0 +1,132 @@
+#!/usr/bin/env python
+import vtk
+from vtk.test import Testing
+from vtk.util.misc import vtkGetDataRoot
+VTK_DATA_ROOT = vtkGetDataRoot()
+
+# The resolution of the density function volume
+res = 100
+
+# Parameters for debugging
+NPts = 1000000
+math = vtk.vtkMath()
+math.RandomSeed(31415)
+
+# create pipeline
+#
+points = vtk.vtkBoundedPointSource()
+points.SetNumberOfPoints(NPts)
+points.ProduceRandomScalarsOn()
+points.ProduceCellOutputOff()
+points.Update()
+
+# Create a sphere implicit function
+sphere = vtk.vtkSphere()
+sphere.SetCenter(0.0,0.1,0.2)
+sphere.SetRadius(0.75)
+
+# Extract points within sphere
+extract = vtk.vtkFitImplicitFunction()
+extract.SetInputConnection(points.GetOutputPort())
+extract.SetImplicitFunction(sphere)
+extract.SetThreshold(0.005)
+extract.GenerateVerticesOn()
+
+# Clip out some of the points with a plane; requires vertices
+plane = vtk.vtkPlane()
+plane.SetOrigin(sphere.GetCenter())
+plane.SetNormal(1,1,1)
+
+clipper = vtk.vtkClipPolyData()
+clipper.SetInputConnection(extract.GetOutputPort())
+clipper.SetClipFunction(plane);
+
+# Generate density field from points
+# Use fixed radius
+dens0 = vtk.vtkPointDensityFilter()
+dens0.SetInputConnection(clipper.GetOutputPort())
+dens0.SetSampleDimensions(res,res,res)
+dens0.SetDensityEstimateToFixedRadius()
+dens0.SetRadius(0.05)
+dens0.SetDensityFormToVolumeNormalized()
+dens0.Update()
+vrange = dens0.GetOutput().GetScalarRange()
+
+map0 = vtk.vtkImageSliceMapper()
+map0.BorderOn()
+map0.SliceAtFocalPointOn()
+map0.SliceFacesCameraOn()
+map0.SetInputConnection(dens0.GetOutputPort())
+
+slice0 = vtk.vtkImageSlice()
+slice0.SetMapper(map0)
+slice0.GetProperty().SetColorWindow(vrange[1]-vrange[0])
+slice0.GetProperty().SetColorLevel(0.5*(vrange[0]+vrange[1]))
+
+# Now densify the point cloud and reprocess
+# Use relative radius
+print("Number of input points: {0}".format(clipper.GetOutput().GetNumberOfPoints()))
+denser = vtk.vtkDensifyPointCloudFilter()
+denser.SetInputConnection(clipper.GetOutputPort())
+denser.SetTargetDistance(0.025)
+denser.SetMaximumNumberOfIterations(5)
+denser.Update()
+print("Number of output points: {0}".format(denser.GetOutput().GetNumberOfPoints()))
+
+dens1 = vtk.vtkPointDensityFilter()
+dens1.SetInputConnection(denser.GetOutputPort())
+dens1.SetSampleDimensions(res,res,res)
+dens1.SetDensityEstimateToFixedRadius()
+dens1.SetRadius(0.05)
+dens1.SetDensityFormToVolumeNormalized()
+dens1.Update()
+vrange = dens1.GetOutput().GetScalarRange()
+
+map1 = vtk.vtkImageSliceMapper()
+map1.BorderOn()
+map1.SliceAtFocalPointOn()
+map1.SliceFacesCameraOn()
+map1.SetInputConnection(dens1.GetOutputPort())
+
+slice1 = vtk.vtkImageSlice()
+slice1.SetMapper(map1)
+slice1.GetProperty().SetColorWindow(vrange[1]-vrange[0])
+slice1.GetProperty().SetColorLevel(0.5*(vrange[0]+vrange[1]))
+
+# Create the RenderWindow, Renderer and both Actors
+#
+ren0 = vtk.vtkRenderer()
+ren0.SetViewport(0,0,0.5,1.0)
+ren1 = vtk.vtkRenderer()
+ren1.SetViewport(0.5,0,1,1.0)
+
+renWin = vtk.vtkRenderWindow()
+renWin.AddRenderer(ren0)
+renWin.AddRenderer(ren1)
+iren = vtk.vtkRenderWindowInteractor()
+iren.SetRenderWindow(renWin)
+
+# Add the actors to the renderer, set the background and size
+#
+ren0.AddActor(slice0)
+ren0.SetBackground(0,0,0)
+ren1.AddActor(slice1)
+ren1.SetBackground(0,0,0)
+
+renWin.SetSize(600,300)
+
+cam = ren0.GetActiveCamera()
+cam.ParallelProjectionOn()
+cam.SetFocalPoint(0,0,0)
+cam.SetPosition(0,0,1)
+ren0.ResetCamera()
+
+ren1.SetActiveCamera(cam)
+
+iren.Initialize()
+
+# render the image
+#
+renWin.Render()
+
+#iren.Start()

Filters/Points/Testing/Python/TestDensifyPointCloudFilter2.py

@@ -0,0 +1,123 @@
+#!/usr/bin/env python
+import vtk
+from vtk.test import Testing
+from vtk.util.misc import vtkGetDataRoot
+VTK_DATA_ROOT = vtkGetDataRoot()
+
+# The resolution of the density function volume
+res = 100
+
+# Parameters for debugging
+NPts = 1000000
+math = vtk.vtkMath()
+math.RandomSeed(31415)
+
+# create pipeline
+#
+# Eight points forming a cube (the regular vtkCubeSource
+# produces duplicate points)
+cube = vtk.vtkPlatonicSolidSource()
+cube.SetSolidTypeToCube()
+
+sphere = vtk.vtkSphereSource()
+sphere.SetRadius(0.05)
+
+glyphs = vtk.vtkGlyph3D()
+glyphs.SetInputConnection(cube.GetOutputPort())
+glyphs.SetSourceConnection(sphere.GetOutputPort())
+glyphs.ScalingOff()
+glyphs.OrientOff()
+
+mapper = vtk.vtkPolyDataMapper()
+mapper.SetInputConnection(glyphs.GetOutputPort())
+
+actor = vtk.vtkActor()
+actor.SetMapper(mapper)
+
+# Now densify the point cloud and reprocess. Use a
+# neighborhood type of N closest points
+denser1 = vtk.vtkDensifyPointCloudFilter()
+denser1.SetInputConnection(cube.GetOutputPort())
+denser1.SetNeighborhoodTypeToNClosest()
+denser1.SetNumberOfClosestPoints(3)
+denser1.SetTargetDistance(1.0)
+denser1.SetMaximumNumberOfIterations(3)
+
+glyphs1 = vtk.vtkGlyph3D()
+glyphs1.SetInputConnection(denser1.GetOutputPort())
+glyphs1.SetSourceConnection(sphere.GetOutputPort())
+glyphs1.ScalingOff()
+glyphs1.OrientOff()
+
+mapper1 = vtk.vtkPolyDataMapper()
+mapper1.SetInputConnection(glyphs1.GetOutputPort())
+
+actor1 = vtk.vtkActor()
+actor1.SetMapper(mapper1)
+
+# Now densify the point cloud and reprocess. Use a
+# neighborhood type of RADIUS
+denser2 = vtk.vtkDensifyPointCloudFilter()
+denser2.SetInputConnection(cube.GetOutputPort())
+denser2.SetNeighborhoodTypeToRadius()
+denser2.SetRadius(1.8)
+denser2.SetTargetDistance(1.0)
+denser2.SetMaximumNumberOfIterations(10)
+denser2.SetMaximumNumberOfPoints(50)
+denser2.Update()
+print(denser2)
+
+glyphs2 = vtk.vtkGlyph3D()
+glyphs2.SetInputConnection(denser2.GetOutputPort())
+glyphs2.SetSourceConnection(sphere.GetOutputPort())
+glyphs2.ScalingOff()
+glyphs2.OrientOff()
+
+mapper2 = vtk.vtkPolyDataMapper()
+mapper2.SetInputConnection(glyphs2.GetOutputPort())
+
+actor2 = vtk.vtkActor()
+actor2.SetMapper(mapper2)
+
+# Create the RenderWindow, Renderer and both Actors
+#
+ren0 = vtk.vtkRenderer()
+ren0.SetViewport(0,0,0.333,1.0)
+ren1 = vtk.vtkRenderer()
+ren1.SetViewport(0.333,0,.6667,1.0)
+ren2 = vtk.vtkRenderer()
+ren2.SetViewport(0.6667,0,1,1.0)
+
+renWin = vtk.vtkRenderWindow()
+renWin.AddRenderer(ren0)
+renWin.AddRenderer(ren1)
+renWin.AddRenderer(ren2)
+iren = vtk.vtkRenderWindowInteractor()
+iren.SetRenderWindow(renWin)
+
+# Add the actors to the renderer, set the background and size
+#
+ren0.AddActor(actor)
+ren0.SetBackground(0,0,0)
+ren1.AddActor(actor1)
+ren1.SetBackground(0,0,0)
+ren2.AddActor(actor2)
+ren2.SetBackground(0,0,0)
+
+renWin.SetSize(900,300)
+
+cam = ren0.GetActiveCamera()
+cam.SetFocalPoint(0,0,0)
+cam.SetPosition(1,1,1)
+ren0.ResetCamera()
+
+ren1.SetActiveCamera(cam)
+ren2.SetActiveCamera(cam)
+
+iren.Initialize()
+
+# render the image
+#
+renWin.Render()
+
+#iren.Start()