Some simulations codes do not operate with grid-based data, but rather with simple particles in 3D space. For example, the Smoothed- Particle Hydrodynamics (SPH) method is a well known method used in fluid dynamics. Besides a very large and better known set of algorithms to visualize the more common grid-based data, ParaView also offers a number of options to deal with particle-based data output. This article presents a review of the techniques we find the most useful, along with the ParaView Python code that would support such visualizations.

Data format issues

ParaView offers native support for some well-known data formats specific to particles. For example, the H5part, AMReX, ENZO, openPMD. Other file formats used in astrophysics for example, Tipsy, Gadget-HDF5, are provided as plugins, developed and maintained on CSCS machines. In some cases, ParaView can also read time-series of such files.

By default, ParaView needs cells, in order to display points. Yet, since particle data are usually not associated with cells, it is often the case that readers will offer the choice of generating such cells. There are two common solutions. A single ParaView cell of type POLYVERTEX can reference millions of particle points, or a cell of type VERTEX can be associated with each point.

Particle queries and selections

As the number of particles grow beyond millions, scientists often want to threshold data based on numerical values, or select sub-sets of particles based on geometric features. We offer here a few examples, extracting particle from a full model:

Extract all particles with rho > 1e-9

Given a ParaView input proxy object called "reader", the following selection is created, and extracted:

selection_sources1 = CreateSelection(proxyname='SelectionQuerySource', registrationname='selection_sources.1', ElementType='Point Data', QueryString='(rho > 1e-9)',
    Selectors=['/Root/Gas'])

selection_filter1 = CreateSelection(proxyname='AppendSelections', registrationname='selection_filter.1', Input=selection_sources1, Expression='s0', SelectionNames=['s0'])

reader = TipsySeriesReader(FileNames=['/local/data/Tipsy/fname'])
reader.PointArrayStatus = ['hsmooth', 'mass', 'rho', 'temp', 'vel']

# create a new 'Extract Selection'
extractSelection1 = ExtractSelection(Input=reader, Selection=selection_filter1)


Extract all particles with rho > 1e-9 AND a negative X coordinate

Qs22 = "(rho > 1e-9) & (points[:,0] < 0)" 

# Note the use of the special keyword 'points' which refers to a numpy-like array of coordinates (x,y, z) tuples.

sel2 = QuerySelect(QueryString=Qs22, Source=reader)

extractSelection2 = ExtractSelection(registrationName='ExtractSelection2', Input=reader)

Extract all particles within a sphere of a given Radius and Center

Center = [-31.310017, -49.53751, -0.11343168]; Radius = 10

Qs3 = "mag(points - [-31.310017, -49.53751, -0.11343168]) < 10"

# note the numpy-like syntax, to compute the magnitude of the coordinates array minus the Center.

sel3 = QuerySelect(QueryString=Qs3, Source=reader)

extractSelection3 = ExtractSelection(registrationName='ExtractSelection3', Input=reader)

Extract all particles within a sphere of a given Radius and Center, using Python Programmable filters

a first example using a geometric object (a sphere)

from vtkmodules.vtkFiltersPoints import vtkExtractEnclosedPoints
from vtkmodules.vtkFiltersSources import vtkSphereSource

Center = [-31.310017, -49.53751, -0.11343168]
Radius = 10

# Create an enclosing surface
ss = vtkSphereSource()
ss.SetPhiResolution(32)
ss.SetThetaResolution(32)
ss.SetCenter(Center)
ss.SetRadius(Radius)

extract = vtkExtractEnclosedPoints()
extract.SetInputData(inputs[0].VTKObject)
extract.SetSurfaceConnection(ss.GetOutputPort())

extract.Update()

output.ShallowCopy(extract.GetOutput())

a second example using an implicit function (a sphere)

from vtkmodules.vtkFiltersPoints import vtkExtractPoints
from vtkmodules.vtkCommonDataModel import vtkSphere
# vtkSphere computes the implicit function and/or gradient for a sphere.
Center = [-31.310017, -49.53751, -0.11343168]
Radius = 10

# Create a sphere implicit function
sphere = vtkSphere()
sphere.SetCenter(Center)
sphere.SetRadius(Radius)

# Extract points within sphere
extract = vtkExtractPoints()
extract.SetInputData(inputs[0].VTKObject)
extract.SetImplicitFunction(sphere)
extract.ExtractInsideOff()
extract.Update()

output.ShallowCopy(extract.GetOutput())

Extract all particles lying within a tolerance off a slicing plane

from vtkmodules.vtkFiltersPoints import vtkFitImplicitFunction
from vtkmodules.vtkCommonDataModel import vtkPlane
# vtkPlane is  a concrete implementation of the abstract class vtkImplicitFunction 
Origin = [36.55, -613.68, -0.14]
Normal = [0, 0, 1]
Tolerance = 0.1
# Create a plane implicit function
plane = vtkPlane()
plane.SetOrigin(Origin)
plane.SetNormal(Normal)

# Extract points "near" the plane, within a tolerance
extract = vtkFitImplicitFunction()
extract.SetInputData(inputs[0].VTKObject)
extract.SetImplicitFunction(plane)
extract.SetThreshold(Tolerance)
extract.Update()
output.ShallowCopy(extract.GetOutput())

Create a hierarchical binning and extract one level of particles

from vtk import vtkHierarchicalBinningFilter, vtkExtractHierarchicalBins

# Bin the points
hBin = vtkHierarchicalBinningFilter()
hBin.SetInputData(inputs[0].VTKObject)
hBin.SetNumberOfLevels(3) # default value

extBin = vtkExtractHierarchicalBins()
extBin.SetInputConnection(hBin.GetOutputPort())
extBin.SetBinningFilter(hBin)
extBin.SetLevel(1)
extBin.SetBin(-1)
extBin.Update()

output.ShallowCopy(extBin.GetOutput())

Create a gridded slicing plane through a particle cloud

ParaView offers two Plane Interpolators. One is generic, using sampling kernels such as Voronoi, Gaussian, or Shepard. The other is SPH-specific, and will use the VTK smoothed-particle hydrodynamics interpolation kernels as described by D.J. Price. Note that the SPH interpolators cannot be used in a mixed mode (MPI + multi-threaded)

# create a new 'SPH Plane Interpolator'
sPHPlaneInterpolator1 = SPHPlaneInterpolator(Input=gas_particles,
   Source='Bounded Plane')
sPHPlaneInterpolator1.DensityArray = 'rho'
sPHPlaneInterpolator1.MassArray = 'Mass'
sPHPlaneInterpolator1.ExcludedArrays = ['Mass']
sPHPlaneInterpolator1.Kernel.SpatialStep = 2.0
sPHPlaneInterpolator1.Source.BoundingBox = [-350.0, 350.0, -350.0, 350.0, -1.0, 1.0]
sPHPlaneInterpolator1.Source.Resolution = N
sPHPlaneInterpolator1.Source.Center = [36.55, -613.68, -0.14]

Reduce the number of particles by Quadric Clustering

ParaView offers two different clustering options. The first one is the standard VTK-based 'Quadric Clustering' and is the same filter used to generate level of detail for ParaView.  It uses a structured grid of bins and merges all points contained in each bin. The second one is VTKm-based 'vtkm Level Of Detail', so it will be fully accelerated (TBB or CUDA) depending on the compilation of ParaView. Because it uses VTKm, it requires a conversion from a VTK dataset, to a VTKm dataset (done internally), at the cost of memory. Further, it does not currently support PolyVertex cells. It will require a vtkPolyData object made of individual vertex cells. If necessary, the ParaView filter "Convert to PointCloud" can be used. Note also that the 'vtkm Level Of Detail' filter will attempt to draw triangles, which does not make sense in the case of particle data. Using a 'Points' representation will turn off this un-desired effect.

# standard VT Quadric Clustering

quadricClustering1 = QuadricClustering( Input=gas)

quadricClustering1.NumberofDimensions = [1024]*3

# VTKm Level of Detail

vTKmLevelOfDetail1 = VTKmLevelOfDetail(Input=convertToPointCloud1)
vTKmLevelOfDetail1.NumberofDimensions = [1024]*3


vTKmLevelOfDetail1Display = Show(vTKmLevelOfDetail1)

vTKmLevelOfDetail1Display.Representation = 'Points'