Is it possible to turn off triangles on vtkClipDataSet filter?

(Bane Sullivan) #1

Is it possible to turn off triangles when clipping a dataset via the vtk.vtkClipDataSet filter? The slicing filter (vtk.vtkCutter) has an awesome GenerateTrianglesOff() method… is there an equivalent for the clip filter?

Box clip with vtkTableBasedClipDataSet and sharp edges?
(Bane Sullivan) #2

Any tips on how to not generate triangles or tetrahedrals when clipping would be greatly appreciated… any ideas?

(Bane Sullivan) #3

Here’s an example from PyVista. In this case, the input mesh (ImageData) has voxels and it gets tessellated into tetras:

import pyvista as pv
from pyvista import examples

mesh = examples.load_uniform()
mesh.plot(show_edges=True)

Then we can use the vtkClipDataSet filter via PyVista’s .clip routine:

clipped = mesh.clip(normal='z')

p = pv.Plotter()
p.add_mesh(clipped, show_edges=True)
p.add_mesh(mesh.outline())
p.show()

Why can’t it remain voxels like how the slice filter can preserve the topology:

slc = mesh.slice()
p = pv.Plotter()
p.add_mesh(slc, show_edges=True)
p.add_mesh(mesh.outline())
p.show()

(Bane Sullivan) #5

@lassoan and I were having a discussion in Opacity issue with point picking on surface that makes more sense to have here

I understand that the point resolution doesn’t change but the cell resolution does - which is really important in my domain (subsurface geoscience/geophysics). Here’s an example: If you have a mesh that is full of voxel-like cells with data attributes on the cell data and you triangulate/tetrahedralize that, then there will be many more cells in the output than the original input mesh. In my domain, we use cell size to represent resolution and if you triangulate a cell then the volume of space now has more data than the original dataset which could be misleading and how much data (information) you have (you know) about that volume of space.

Here’s an example mesh: https://www.dropbox.com/s/5ivldekj8weqc79/treemesh.vtk?dl=0

and some example code to demo this issue using PyVista:

import pyvista as pv

mesh = pv.read('/Users/bane/Desktop/treemesh.vtk')
subsurface = mesh.threshold(0.5, scalars='Active')

subsurface.plot(show_edges=True, color='w')

Then if run a vtkDataSetTriangleFilter, I get this:

subsurface.triangulate().plot(show_edges=True, color='w')

or a vtkClipDataSet filter which generates triangles, I get this:

subsurface.clip().plot(show_edges=True, color='w',)

Both of which are misleading to the actual resolution we are able to recover from our physical modeling techniques. So if I wanted to extract a region and share it with a colleage, then they would receive a misrepresentation of the resolution at which I can investigate the subsurface.

Which is why I care so much about this turning off triangles when clipping a dataset - I see no reason why it isn’t possible

(Andras Lasso) #6

Clipping is an extremely complicated operation. It is hard enough to implement it for triangles or tetrahedral elements, so I guess that’s why the mesh is converted to these elements before clipping. Maybe you could give it a try and implement clipping of cell types that are important for you in VTK. You could also add pedigree IDs to cells so that you can trace back which original cell each new cell corresponds to.

From the examples above it is not clear at all what’s happening - what VTK filters are used and how. If you would like to have VTK users or developers to look at a problem, then probably you need to provide VTK examples.

(Bane Sullivan) #7

Fair. There’s a lot of other things going on like meshio and rendering (mappers/actors/etc.) and I figured I’d keep the code simple but here is exactly what is going on for the .clip() filter which leverages the vtkClipDataSet.

import vtk
import numpy as np

NORMALS = {
    'x': [1, 0, 0],
    'y': [0, 1, 0],
    'z': [0, 0, 1],
    '-x': [-1, 0, 0],
    '-y': [0, -1, 0],
    '-z': [0, 0, -1],
}

def clip(dataset, normal='x', origin=None, invert=True):
    """
    Clip a dataset by a plane by specifying the origin and normal. If no
    parameters are given the clip will occur in the center of that dataset

    Parameters
    ----------
    normal : tuple(float) or str
        Length 3 tuple for the normal vector direction. Can also be
        specified as a string conventional direction such as ``'x'`` for
        ``(1,0,0)`` or ``'-x'`` for ``(-1,0,0)``, etc.

    origin : tuple(float)
        The center ``(x,y,z)`` coordinate of the plane on which the clip
        occurs

    invert : bool
        Flag on whether to flip/invert the clip

    """
    if isinstance(normal, str):
        normal = NORMALS[normal.lower()]
    # find center of data if origin not specified
    if origin is None:
        origin = dataset.GetCenter()
    # create the plane for clipping
    plane = vtk.vtkPlane()
    # NORMAL MUST HAVE MAGNITUDE OF 1
    normal = normal / np.linalg.norm(normal)
    plane.SetNormal(normal)
    plane.SetOrigin(origin)
    # run the clip
    if isinstance(dataset, vtk.vtkPolyData):
        alg = vtk.vtkClipPolyData()
    elif isinstance(dataset, vtk.vtkImageData):
        alg = vtk.vtkClipVolume()
        alg.SetMixed3DCellGeneration(True)
    else:
        alg = vtk.vtkClipDataSet()
    alg.SetInputDataObject(dataset) # Use the grid as the data we desire to cut
    alg.SetClipFunction(plane) # the the cutter to use the plane we made
    alg.SetInsideOut(invert) # invert the clip if needed
    alg.Update() # Perfrom the Cut
    return alg.GetOutput()

In the example above, I’m passing a vtkUnstructuredGrid to the function

(Bane Sullivan) #8

That actually might work… perhaps I could keep track of the original cell ID’s and then use those to extract the cells from the original mesh using a selection routine with vtkSelectionNode, vtkSelection, and vtkExtractSelection. That parts easy… but how do I keep track of the original cell IDs when using the vtkClipDataSet filter?

(Bane Sullivan) #9

I could add the cell IDs to the cell data of the input then catch them on the output :man_facepalming: :bulb:

(Bane Sullivan) #10

Only issue is that the clipping plane wouldn’t be preserved in the output and I wouldn’t be able to make intersection polygons:

...
ids = np.arange(0, subsurface.n_cells, dtype=int)
subsurface['ids'] = ids
clp = subsurface.clip([-1,1,1])
clp.plot(show_edges=True, color='w', use_panel=False)

clp_preserved = subsurface.extract_cells(np.unique(clp['ids']))
clp_preserved.plot(show_edges=True, color='w', use_panel=False)

(Bane Sullivan) #11

Hm… maybe I could use the slice filter and leverage the fact that vtk.vtkCutter has an awesome GenerateTrianglesOff() method to make those intersection polygons then merge the slice and the clipped mesh using IDs to make a hack of a solution. But I think it’ll work and could be generalized… I’ll try it!

(Bill Lorensen) #12

ClipPolyData does preserve original polygons that are not clipped. Only triangles if clipped. I guess you probably know that.
See this example.

2 Likes
(Bane Sullivan) #13

Yeah, the ClipPolyData filter is great for 2D cells - but I’m mostly clipping volumetric cells in rectilinear grids, un/structured grids, and image data so having a clear way to clip and preserve topology not intersected by the clip would be huge for me

(Bill Lorensen) #14

Have you tried vtkGenericClip?

(Andras Lasso) #15

When you clip volumetric meshes, you have the option of not cutting into the cells, just keep or remove whole cells.

(Bill Lorensen) #16

I just checked, GenericClip also generates triangles…

1 Like
(Bill Lorensen) #17

Could you share a sample dataset? I may have a solution using vtkMultiThreshold. I can extract the inside and border voxels as separate datasets. hen clip the border cells and append them with the inside.

I’m not done yet, but here is the result of the MultiThreshold.

(Bane Sullivan) #18

@lorensen, that looks awesome! Here is an unstructured grid that I have been working with lately where the cell resolution is really important: https://www.dropbox.com/s/5ivldekj8weqc79/treemesh.vtk?dl=0

It’d be awesome if we could find a way to preserve cells not touched by the clipping plane and generate intersection polygons where the plane intersect cells in the mesh.

(Vincent Rivola) #19

Hi,

Did you try vtkTableBasedClipDataSet ?
From documentation it says:
" It is worth mentioning that vtkTableBasedClipDataSet is capable of addressing the artifacts that may occur with vtkClipDataSet due to the possibly inconsistent triangulation modes between neighboring cells. In addition, the former is much faster than the latter. Furthermore, the former produces less cells (with ratio usually being 5~6) than by the latter in the output. In other words, this filter retains the original cells (i.e., without triangulation / tetrahedralization) wherever possible"

1 Like
(Bill Lorensen) #20

Vincent is correct. This example shows that triangles are NOT produced as shown in this closeup.

. Bane, ill try it with your dataset later today.

1 Like
(Bill Lorensen) #21

And, This example clips a rectilinear grid. It also reports the cell types of the clipped data.

The clipped dataset(inside) contains a 
vtkUnstructuredGrid that has 14176 cells
Cell type vtkTetra occurs 3084 times.
Cell type vtkHexahedron occurs 6820 times.
Cell type vtkWedge occurs 1196 times.
Cell type vtkPyramid occurs 3076 times.

The clipped dataset(outside) contains a 
vtkUnstructuredGrid that has 125600 cells
Cell type vtkTetra occurs 3276 times.
Cell type vtkHexahedron occurs 117732 times.
Cell type vtkWedge occurs 1260 times.
Cell type vtkPyramid occurs 3332 times.

1 Like