Locate cell containing a point


I’m trying to locate a point (with xyz coordinates) in a grid (i.e. find the cell containing that point). It can be any type of grid: image data, structured, explicit structured, unstructured. My goal is to get interpolation weigths to interpolate vertex data.

I’m using the FindCell method of vtkDataSet (https://vtk.org/doc/nightly/html/classvtkDataSet.html#a2221c10d3c4cca44e82c5ef70e4e1cbd, see code below).

It works fine when the dataset is an image data, but for StructuredGrids, I have a crash (see stack trace below).

Am I doing something wrong? Or is there a better way to do what I’m trying to do?

Any help will be greatly appreciated




// Get the cell containing the point
int sub_id;
double pcoords[3];
vtkIdType cell_id = dataSet->FindCell(point, nullptr, -1, 0, sub_id, pcoords, nullptr);

//Get the number of points of the cell
vtkNew<vtkIdList> cell_points;
dataSet->GetCellPoints(cell_id, cell_points);
vtkIdType nb_points = cell_points->GetNumberOfIds();

//Get the interpolation weights
dataSet->FindCell(point, nullptr, cell_id, 0, sub_id, pcoords, weights);

Stacktrace error:
… my method
7 0x7fe7ebe82c61 vtkPointSet::FindCell(double*, vtkCell*, long long, double, int&, double*, double*) + 99
6 0x7fe7ebe82bb2 vtkPointSet::FindCell(double*, vtkCell*, vtkGenericCell*, long long, double, int&, double*, double*) + 188
5 0x7fe7ebc96afd vtkClosestPointStrategy::FindCell(double*, vtkCell*, vtkGenericCell*, long long, double, int&, double*, double*) + 805
4 0x7fe7ebc967af /home/fuet/SPDE/ParaView/ParaView57/build/lib/libvtkCommonDataModel-pv5.7.so.1(+0x1e97af) [0x7fe7ebc967af]
3 0x7fe7ebc96624 /home/fuet/SPDE/ParaView/ParaView57/build/lib/libvtkCommonDataModel-pv5.7.so.1(+0x1e9624) [0x7fe7ebc96624]
2 0x7fe7ebd13a90 vtkHexahedron::EvaluatePosition(double const*, double*, int&, double*, double&, double*) + 738
1 0x7fe7ebd14382 vtkHexahedron::InterpolationFunctions(double const*, double*) + 114
0 0x7fe7f4153f20 /lib/x86_64-linux-gnu/libc.so.6(+0x3ef20) [0x7fe7f4153f20]
( 29.308s) [paraview ] :0 FATL| Signal: SIGSEGV
Erreur de segmentation (core dumped)

The answer to this question is here: