Probe an unstructured grid (vtkUnstructuredGrid) with a volume (three-dimensional vtkImageData)

Hi everyone!

I’m trying to use the vtkProbeFilter to interpolate some PointData on a mesh to the points of an Image. As specified in the doc, I believe it is possible:

For example, an unstructured grid (vtkUnstructuredGrid) can be probed with a volume (three-dimensional vtkImageData)

However, I’m wondering if I’m doing it the right way ; I’m first converting my unstructured grid to a polydata because it seems to be the way to further create the image (input of the filter). For example I didn’t find a vtkUnstructuredGridToImageStencil class…

But having my ImageData (which is valid), my problem is in the filter. Its output is empty, the file only contains the header, see file below.

Do you have any clue of why it’s empty? Thanks a lot to take the time to answer me.
I’m working with the last version of VTK (vtk-9.0.0).

Here is my code (I omitted the UnstructuredGridReader and Image conversion):

vtkSmartPointer<vtkUnstructuredGrid> source; // Previously read from file, refered as the mesh
// Convert UnstructuredGrid to Polydata with GeometryFilter
// Convert PolyData to ImageData with PolyDataToImageStencil and ImageStencil
vtkSmartPointer<vtkImageData> input; // Image of the mesh

// Add PointData to the mesh
vtkSmartPointer<vtkFloatArray> ptData = vtkSmartPointer<vtkFloatArray>::New();
for(Int i=0; i<source->GetNumberOfPoints(); i++) ptData->InsertTuple1(i, foo);
source->GetPointData()->AddArray(ptData);

// Interpolate PointData from the mesh to the points of the Image
vtkSmartPointer<vtkProbeFilter> filter = vtkSmartPointer<vtkProbeFilter>::New();
filter->SetSourceData(source);
filter->SetInputData(input);
filter->Update();

// Write image to file
vtkSmartPointer<vtkMetaImageWriter> writer = vtkSmartPointer<vtkMetaImageWriter>::New();
writer->SetFileName("image.mha");
writer->SetInputData(filter->GetOutput());
writer->Write();

And my image.mha file:

ObjectType = Image
NDims = 3
BinaryData = True
BinaryDataByteOrderMSB = False
CompressedData = True
TransformMatrix = 1 0 0 0 1 0 0 0 1
Offset = 0 0 0
CenterOfRotation = 0 0 0
ElementSpacing = 1 1 1
DimSize = 124 112 110
AnatomicalOrientation = ???
ElementType = MET_DOUBLE
ElementDataFile = LOCAL

You have forgot to set the array name. Unnamed arrays cannot be added to a data set.

To VTK developers: New VTK users often make the mistake of adding an array without specifying name first. In general, you don’t need to set default values in VTK, which makes VTK easy to use. This is an exception - where you need to set a default name manually, even if you don’t care what that name is. Would it be possible to set a default array name when an array is created or when an unnamed array is added to a data set?

Indeed, there is no reason not to name an array by default. An empty string or “Array” should do the trick.

That’s not strictly true. You can add an array without a name. It’s just not particularly useful. You can indeed set an unnamed array as an active attribute (e.g., vtkPointData::SetScalars(), vtkPointData::SetVectors(), etc.) and filters that use the attribute API can use them, but that’s old-style VTK and nowadays naming the array to use is strongly preferred.

One problem is that array names are required to be unique (see the implementation of vtkFieldData::SetArray(int i, vtkAbstractArray* data), so if you added a default, you would need to make it unique from names of other arrays already present.

At many places, VTK assumes that array name is set and otherwise behaves unexpected ways or it may even crash.

For example, I remember that vtkTable can crash if it contains unnamed array, and there are a number of small things, such as you can add an unnamed array many times to a vtkPointData (the mechanism that would prevent this relies on array name), vtkArrayData::GetArrayByName crashes if any of the arrays are unnamed, vtkGenericAttribute::PrintSelf crashes, etc.

If the decision is that unnamed arrays are supported then these issues and inconsistencies would need to be fixed (should not be difficult, there could be just a number of them hiding here and there). If unnamed arrays are officially not supported or deprecated then there should be at least warning messages for example when they are added to a data set.

I agree, we should require arrays to be named and filters should all operate by array names. Active attributes should be completely phased out. We are slowly moving that way I think, but it is a big job to go through all of VTK and make that happen.

1 Like

I’m also using VTKProbeFilter to do the same thing and try to write the output as a .mhd file, but the output file is empty and only contains the header, see below.
Here is my code and can you please let me know any issues regarding the implementation?

Am I doing something wrong here ? Please help me on this matter.

vtkSmartPointer<vtkUnstructuredGridReader> reader =
        vtkSmartPointer<vtkUnstructuredGridReader>::New();
reader->SetFileName(ref_volumetric_mesh_path.c_str());
reader->ReadAllVectorsOn();
reader->ReadAllScalarsOn();
reader->Update();
reader->GetOutput()->GetScalarRange();

vtkSmartPointer<vtkUnstructuredGrid> source = reader->GetOutput();

vtkSmartPointer<vtkImageData> input =
	vtkSmartPointer<vtkImageData>::New();

input->SetDimensions(60, 60, 100);
input->SetSpacing(0.5, 0.5, 0.5);
input->SetOrigin(0.0, 0.0, 0.0);

input->AllocateScalars(VTK_DOUBLE, 1);

// Interpolate PointData from the mesh to the points of the Image
vtkSmartPointer<vtkProbeFilter> filter = vtkSmartPointer<vtkProbeFilter>::New();
filter->SetInputData(input);
filter->SetSourceData(source);
filter->Update();

vtkSmartPointer<vtkMetaImageWriter> writer = vtkSmartPointer<vtkMetaImageWriter>::New();
writer->SetFileName(ref_volumetric_mesh_voxelized_path.c_str());
writer->SetInputData(filter->GetOutput());
writer->Write();

.mhd file only contains the header. See below.

    ObjectType = Image
    NDims = 3
    BinaryData = True
    BinaryDataByteOrderMSB = False
    CompressedData = True
    TransformMatrix = 1 0 0 0 1 0 0 0 1
    Offset = 0 0 0
    CenterOfRotation = 0 0 0
    ElementSpacing = 0.5 0.5 0.5
    DimSize = 60 60 100
    AnatomicalOrientation = ???
    ElementType = MET_DOUBLE
    ElementDataFile = example_gen_image_vol.zraw
1 Like

Hi, I wonder did you get this problem solved?
I am having touble to use VTKprobeFilter with unstructured grid data.
if you have any suggestion, could you please share your thought?