[vtkusers] vtkConnectivityFilter on vtkStructuredGrid

[David Doria]

Hi all,

The ConnectivityFilter seems to work as expected with PolyData (e.g.
http://www.vtk.org/Wiki/VTK/Examples/Cxx/Filtering/ConnectivityFilter).
However, I've been trying to find "connected" regions (obviously by
some property, as the whole grid is topologically connected) in a
vtkStructuredGrid and it does not seem to respect the specified
ScalarConnectivity. As a minimal example, I created a grid with the
top and bottom slices set to "valid" (an arbitrarily chosen value of
10) and the interior slices set to "invalid" (0). I used my own
vtkIntArray instead of the vtkGhostType set by blanking just to be
sure it wasn't any weirdness with that. I set this array as the cell
scalars of the grid, specified to use it with:

connectivityFilter->SetInputArrayToProcess(0, 0, 0,
vtkDataObject::FIELD_ASSOCIATION_CELLS, "ValidValue");

and then used:

    connectivityFilter->ScalarConnectivityOn();
    connectivityFilter->SetScalarRange(9,11);

to try to find regions where the value is 10. In this case, I'd expect
two regions (one for each layer that is "valid"). However, the filter
outputs only one region (the entire grid). The sample code is below -
can anyone see what is going wrong?

#include <vtkCellData.h>
#include <vtkIntArray.h>
#include <vtkConnectivityFilter.h>
#include <vtkDataSetAttributes.h>
#include <vtkSmartPointer.h>
#include <vtkStructuredGrid.h>

int main(int, char *[])
{
    // Create a 7x7 grid, with only 2 unblanked z slices
    vtkSmartPointer<vtkStructuredGrid> structuredGrid =
      vtkSmartPointer<vtkStructuredGrid>::New();

    vtkSmartPointer<vtkPoints> points =
      vtkSmartPointer<vtkPoints>::New();

    unsigned int voxelsPerDimension = 7;

    // Create an array
    vtkSmartPointer<vtkIntArray> valid =
      vtkSmartPointer<vtkIntArray>::New();
    valid->SetName("ValidValue");
    valid->SetNumberOfComponents(1);
    valid->SetNumberOfTuples(voxelsPerDimension * voxelsPerDimension *
voxelsPerDimension);

    // The grid is created as points (the corners of voxels), so we
need n+1 points per dimension
    for(unsigned int k = 0; k < voxelsPerDimension + 1; k++)
    {
        for(unsigned int j = 0; j < voxelsPerDimension + 1; j++)
        {
          for(unsigned int i = 0; i < voxelsPerDimension + 1; i++)
          {
              points->InsertNextPoint(i, j, k);
          }
        }
    }

    // Specify the dimensions of the grid (again, using point dimensions)
    structuredGrid->SetDimensions(voxelsPerDimension + 1,
voxelsPerDimension + 1, voxelsPerDimension + 1);
    structuredGrid->SetPoints(points);

    int dims[3];
    structuredGrid->GetDimensions(dims);

    for(unsigned int i = 0; i < voxelsPerDimension; ++i)
    {
        for(unsigned int j = 0; j < voxelsPerDimension; ++j)
        {
            for(unsigned int k = 0; k < voxelsPerDimension; ++k)
            {
                int ijk[3];
                ijk[0] = i;
                ijk[1] = j;
                ijk[2] = k;

                vtkIdType cellID = vtkStructuredData::ComputeCellId(dims,ijk);

                if(k == 0 || k == 6)
                {
                    valid->SetTuple1(cellID, 10);
                }
                else
                {
                    structuredGrid->BlankCell(cellID);
                    valid->SetTuple1(cellID, 0);
                }

            }
        }
    }

    structuredGrid->GetCellData()->SetScalars(valid);
    structuredGrid->Modified();

    vtkSmartPointer<vtkConnectivityFilter> connectivityFilter =
vtkConnectivityFilter::New();
    connectivityFilter->SetInputArrayToProcess(0, 0, 0,
vtkDataObject::FIELD_ASSOCIATION_CELLS, "ValidValue");
    connectivityFilter->SetInputData(structuredGrid);
    connectivityFilter->SetExtractionModeToAllRegions();
    connectivityFilter->ColorRegionsOn();
    connectivityFilter->ScalarConnectivityOn();
    connectivityFilter->SetScalarRange(9,11);
    connectivityFilter->Update();

    std::cout << "There were " <<
connectivityFilter->GetNumberOfExtractedRegions() << " extracted
regions." << std::endl;

    return EXIT_SUCCESS;
}

Thanks,

David