[Resolved] How to create a template from ConvertToVTK example?

Hello everyone,

I work on Mac OS X 10.13.6 with Clang 9.1 (Apple latest version) and ITK v 4.13.0.

I try to “templatize” the example ConvertToVTK.

To do this, I just comment the first typedef and add some template definitions to VisitVTKCellsClass class, ConvertMeshToUnstructuredGrid and CreateMeshWithEdges functions.

It’s the function ConvertMeshToUnstructuredGrid that pose problem: I can’t add the visitors tv, qv and lv because of some conversions problems.

The following code is my failed attempt…:

Summary

// ITK

#include <itkLineCell.h>

#include <itkMesh.h>

#include <itkTriangleCell.h>

#include <itkQuadrilateralCell.h>

// VTK

#include “vtkVersion.h”

#include <vtkCellArray.h>

#include <vtkSmartPointer.h>

#include <vtkUnstructuredGrid.h>

#include <vtkXMLUnstructuredGridWriter.h>

// Functions

template <typename MeshType>

static typename MeshType::Pointer CreateMeshWithEdges();

template <typename MeshType>

static void ConvertMeshToUnstructuredGrid(typename MeshType::Pointer, vtkUnstructuredGrid *);

template<typename MeshType>

class VisitVTKCellsClass

{

vtkCellArray *m_Cells;

int *m_LastCell;

int *m_TypeArray;

public:

// typedef the itk cells we are interested in

typedef itk::CellInterface<

typename MeshType::PixelType,

typename MeshType::CellTraits>

CellInterfaceType;

typedef itk::LineCell<CellInterfaceType> floatLineCell;

typedef itk::TriangleCell<CellInterfaceType> floatTriangleCell;

typedef itk::QuadrilateralCell<CellInterfaceType> floatQuadrilateralCell;

// Set the vtkCellArray that will be constructed

void SetCellArray(vtkCellArray *a)

{

m_Cells = a;

}

// Set the cell counter pointer

void SetCellCounter(int *i)

{

m_LastCell = i;

}

// Set the type array for storing the vtk cell types

void SetTypeArray(int *i)

{

m_TypeArray = i;

}

// Visit a triangle and create the VTK_TRIANGLE cell

void Visit(unsigned long, floatTriangleCell *t)

{

m_Cells->InsertNextCell(3, (vtkIdType *)t->PointIdsBegin());

m_TypeArray[*m_LastCell] = VTK_TRIANGLE;

(*m_LastCell)++;

}

// Visit a triangle and create the VTK_QUAD cell

void Visit(unsigned long, floatQuadrilateralCell *t)

{

m_Cells->InsertNextCell(4, (vtkIdType *)t->PointIdsBegin());

m_TypeArray[*m_LastCell] = VTK_QUAD;

(*m_LastCell)++;

}

// Visit a line and create the VTK_LINE cell

void Visit(unsigned long, floatLineCell *t)

{

m_Cells->InsertNextCell(2, (vtkIdType *)t->PointIdsBegin());

m_TypeArray[*m_LastCell] = VTK_LINE;

(*m_LastCell)++;

}

};

int main(int, char *[])

{

typedef itk::Mesh<double, 3> MeshType;

MeshType::Pointer mesh = CreateMeshWithEdges<MeshType>();

vtkSmartPointer<vtkUnstructuredGrid> unstructuredGrid = vtkSmartPointer<vtkUnstructuredGrid>::New();

ConvertMeshToUnstructuredGrid<MeshType>(mesh, unstructuredGrid);

// Write file

vtkSmartPointer<vtkXMLUnstructuredGridWriter> writer =

vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();

writer->SetFileName(“output.vtu”);

#if VTK_MAJOR_VERSION <= 5

writer->SetInputConnection(unstructuredGrid->GetProducerPort());

#else

writer->SetInputData(unstructuredGrid);

#endif

writer->Write();

return EXIT_SUCCESS;

}

template <typename MeshType>

typename MeshType::Pointer CreateMeshWithEdges()

{

typename MeshType::Pointer mesh = MeshType::New();

// Create 4 points and add them to the mesh

typename MeshType::PointType p0, p1, p2, p3;

// the same code from the example

}

template <typename MeshType>

void ConvertMeshToUnstructuredGrid(typename MeshType::Pointer mesh, vtkUnstructuredGrid *unstructuredGrid)

{

// Get the number of points in the mesh

int numPoints = mesh->GetNumberOfPoints();

if (numPoints == 0)

{

mesh->Print(std::cerr);

std::cerr << "no points in Grid " << std::endl;

exit(-1);

}

// Create the vtkPoints object and set the number of points

vtkPoints *vpoints = vtkPoints::New();

vpoints->SetNumberOfPoints(numPoints);

// Iterate over all the points in the itk mesh filling in

// the vtkPoints object as we go

typename MeshType::PointsContainer::Pointer points = mesh->GetPoints();

// In ITK the point container is not necessarily a vector, but in VTK it is

vtkIdType VTKId = 0;

std::map<vtkIdType, int> IndexMap;

for (typename MeshType::PointsContainer::Iterator i = points->Begin();

i != points->End(); ++i, VTKId++)

{

// Get the point index from the point container iterator

IndexMap[VTKId] = i->Index();

// Set the vtk point at the index with the the coord array from itk

// itk returns a const pointer, but vtk is not const correct, so

// we have to use a const cast to get rid of the const

vpoints->SetPoint(VTKId, const_cast<float *>(i->Value().GetDataPointer()));

}

// Set the points on the vtk grid

unstructuredGrid->SetPoints(vpoints);

// Setup some VTK things

int vtkCellCount = 0; // running counter for current cell being inserted into vtk

int numCells = mesh->GetNumberOfCells();

int *types = new int[numCells]; // type array for vtk

// create vtk cells and estimate the size

vtkCellArray *cells = vtkCellArray::New();

cells->EstimateSize(numCells, 4);

// Setup the line visitor

typedef itk::CellInterfaceVisitorImplementation<

float, typename MeshType::CellTraits,

itk::LineCell<itk::CellInterface<typename MeshType::PixelType, typename MeshType::CellTraits>>,

VisitVTKCellsClass<MeshType>>

LineVisitor;

typename LineVisitor::Pointer lv = LineVisitor::New();

lv->SetTypeArray(types);

lv->SetCellCounter(&vtkCellCount);

lv->SetCellArray(cells);

// Setup the triangle visitor

typedef itk::CellInterfaceVisitorImplementation<

float, typename MeshType::CellTraits,

itk::TriangleCell<itk::CellInterface<typename MeshType::PixelType, typename MeshType::CellTraits>>,

VisitVTKCellsClass<MeshType>>

TriangleVisitor;

typename TriangleVisitor::Pointer tv = TriangleVisitor::New();

tv->SetTypeArray(types);

tv->SetCellCounter(&vtkCellCount);

tv->SetCellArray(cells);

// Setup the quadrilateral visitor

typedef itk::CellInterfaceVisitorImplementation<

float, typename MeshType::CellTraits,

itk::QuadrilateralCell<itk::CellInterface<typename MeshType::PixelType, typename MeshType::CellTraits>>,

VisitVTKCellsClass<MeshType>>

QuadrilateralVisitor;

typename QuadrilateralVisitor::Pointer qv = QuadrilateralVisitor::New();

qv->SetTypeArray(types);

qv->SetCellCounter(&vtkCellCount);

qv->SetCellArray(cells);

// Add the visitors to a multivisitor

typename MeshType::CellType::MultiVisitor::Pointer mv = MeshType::CellType::MultiVisitor::New();

// error: no viable conversion from…
mv->AddVisitor(tv);

mv->AddVisitor(qv);

mv->AddVisitor(lv);

// Now ask the mesh to accept the multivisitor which

// will Call Visit for each cell in the mesh that matches the

// cell types of the visitors added to the MultiVisitor

mesh->Accept(mv);

// Now set the cells on the vtk grid with the type array and cell array

unstructuredGrid->SetCells(types, cells);

std::cout << “Unstructured grid has " << unstructuredGrid->GetNumberOfCells() << " cells.” << std::endl;

// Clean up vtk objects

cells->Delete();

vpoints->Delete();

}

You defined the mesh type with double, and you were using floats everywhere else. Here is a modified code which does compile for me.

// ITK
#include <itkLineCell.h>
#include <itkMesh.h>
#include <itkTriangleCell.h>
#include <itkQuadrilateralCell.h>
// VTK
#include <vtkVersion.h>
#include <vtkCellArray.h>
#include <vtkSmartPointer.h>
#include <vtkUnstructuredGrid.h>
#include <vtkXMLUnstructuredGridWriter.h>

template<typename MeshType>
class VisitVTKCellsClass
{
    vtkCellArray *m_Cells;
    int *m_LastCell;
    int *m_TypeArray;
public:
    // typedef the itk cells we are interested in
    typedef itk::CellInterface<
        typename MeshType::PixelType,
        typename MeshType::CellTraits>
        CellInterfaceType;
    typedef itk::LineCell<CellInterfaceType> floatLineCell;
    typedef itk::TriangleCell<CellInterfaceType> floatTriangleCell;
    typedef itk::QuadrilateralCell<CellInterfaceType> floatQuadrilateralCell;
    // Set the vtkCellArray that will be constructed
    void SetCellArray(vtkCellArray *a)
    {
        m_Cells = a;
    }
    // Set the cell counter pointer
    void SetCellCounter(int *i)
    {
        m_LastCell = i;
    }
    // Set the type array for storing the vtk cell types
    void SetTypeArray(int *i)
    {
        m_TypeArray = i;
    }
    // Visit a triangle and create the VTK_TRIANGLE cell
    void Visit(unsigned long, floatTriangleCell *t)
    {
        m_Cells->InsertNextCell(3, (vtkIdType *)t->PointIdsBegin());
        m_TypeArray[*m_LastCell] = VTK_TRIANGLE;
        (*m_LastCell)++;
    }
    // Visit a triangle and create the VTK_QUAD cell
    void Visit(unsigned long, floatQuadrilateralCell *t)
    {
        m_Cells->InsertNextCell(4, (vtkIdType *)t->PointIdsBegin());
        m_TypeArray[*m_LastCell] = VTK_QUAD;
        (*m_LastCell)++;
    }
    // Visit a line and create the VTK_LINE cell
    void Visit(unsigned long, floatLineCell *t)
    {
        m_Cells->InsertNextCell(2, (vtkIdType *)t->PointIdsBegin());
        m_TypeArray[*m_LastCell] = VTK_LINE;
        (*m_LastCell)++;
    }
};


template <typename MeshType>
typename MeshType::Pointer CreateMeshWithEdges()
{
    typename MeshType::Pointer mesh = MeshType::New();
    // Create 4 points and add them to the mesh
    typename MeshType::PointType p0, p1, p2, p3;
    // the same code from the example
    return mesh;
}

template <typename MeshType>
void ConvertMeshToUnstructuredGrid(typename MeshType::Pointer mesh, vtkUnstructuredGrid *unstructuredGrid)
{
    // Get the number of points in the mesh
    int numPoints = mesh->GetNumberOfPoints();
    if (numPoints == 0)
    {
        mesh->Print(std::cerr);
        std::cerr << "no points in Grid " << std::endl;
        exit(-1);
    }
    // Create the vtkPoints object and set the number of points
    vtkPoints *vpoints = vtkPoints::New();
    vpoints->SetNumberOfPoints(numPoints);
    // Iterate over all the points in the itk mesh filling in
    // the vtkPoints object as we go
    typename MeshType::PointsContainer::Pointer points = mesh->GetPoints();
    // In ITK the point container is not necessarily a vector, but in VTK it is
    vtkIdType VTKId = 0;
    std::map<vtkIdType, int> IndexMap;
    for (typename MeshType::PointsContainer::Iterator i = points->Begin();
        i != points->End(); ++i, VTKId++)
    {
        // Get the point index from the point container iterator
        IndexMap[VTKId] = i->Index();
        // Set the vtk point at the index with the the coord array from itk
        // itk returns a const pointer, but vtk is not const correct, so
        // we have to use a const cast to get rid of the const
        vpoints->SetPoint(VTKId, const_cast<float *>(i->Value().GetDataPointer()));
    }
    // Set the points on the vtk grid
    unstructuredGrid->SetPoints(vpoints);
    // Setup some VTK things
    int vtkCellCount = 0; // running counter for current cell being inserted into vtk
    int numCells = mesh->GetNumberOfCells();
    int *types = new int[numCells]; // type array for vtk
                                    // create vtk cells and estimate the size
    vtkCellArray *cells = vtkCellArray::New();
    cells->EstimateSize(numCells, 4);
    // Setup the line visitor
    typedef itk::CellInterfaceVisitorImplementation<
        typename MeshType::PixelType, typename MeshType::CellTraits,
        itk::LineCell<typename MeshType::CellType>,
        VisitVTKCellsClass<MeshType>>
        LineVisitor;
    typename LineVisitor::Pointer lv = LineVisitor::New();
    lv->SetTypeArray(types);
    lv->SetCellCounter(&vtkCellCount);
    lv->SetCellArray(cells);
    // Setup the triangle visitor
    typedef itk::CellInterfaceVisitorImplementation<
        typename MeshType::PixelType, typename MeshType::CellTraits,
        itk::TriangleCell<typename MeshType::CellType>,
        VisitVTKCellsClass<MeshType>>
        TriangleVisitor;
    typename TriangleVisitor::Pointer tv = TriangleVisitor::New();
    tv->SetTypeArray(types);
    tv->SetCellCounter(&vtkCellCount);
    tv->SetCellArray(cells);
    // Setup the quadrilateral visitor
    typedef itk::CellInterfaceVisitorImplementation<
        typename MeshType::PixelType, typename MeshType::CellTraits,
        itk::QuadrilateralCell<typename MeshType::CellType>,
        VisitVTKCellsClass<MeshType>>
        QuadrilateralVisitor;
    typename QuadrilateralVisitor::Pointer qv = QuadrilateralVisitor::New();
    qv->SetTypeArray(types);
    qv->SetCellCounter(&vtkCellCount);
    qv->SetCellArray(cells);
    // Add the visitors to a multivisitor
    typename MeshType::CellType::MultiVisitor::Pointer mv = MeshType::CellType::MultiVisitor::New();
    // error: no viable conversion from…
    mv->AddVisitor(tv);
    mv->AddVisitor(qv);
    mv->AddVisitor(lv);
    // Now ask the mesh to accept the multivisitor which
    // will Call Visit for each cell in the mesh that matches the
    // cell types of the visitors added to the MultiVisitor
    mesh->Accept(mv);
    // Now set the cells on the vtk grid with the type array and cell array
    unstructuredGrid->SetCells(types, cells);
    std::cout << "Unstructured grid has " << unstructuredGrid->GetNumberOfCells() << " cells." << std::endl;
    // Clean up vtk objects
    cells->Delete();
    vpoints->Delete();
}

int main(int, char *[])
{
    typedef itk::Mesh<double, 3> MeshType;
    MeshType::Pointer mesh = CreateMeshWithEdges<MeshType>();
    vtkSmartPointer<vtkUnstructuredGrid> unstructuredGrid = vtkSmartPointer<vtkUnstructuredGrid>::New();
    ConvertMeshToUnstructuredGrid<MeshType>(mesh, unstructuredGrid);
    // Write file
    vtkSmartPointer<vtkXMLUnstructuredGridWriter> writer =
        vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();
    writer->SetFileName("output.vtu");
#if VTK_MAJOR_VERSION <= 5
    writer->SetInputConnection(unstructuredGrid->GetProducerPort());
#else
    writer->SetInputData(unstructuredGrid);
#endif
    writer->Write();
    return EXIT_SUCCESS;
}

Thank you so much! It works for me too.
I didn’t see the remaining float…

Is it possible to put this generic example on the website instead of the previous one ?

It can be quite useful for many people, don’t you think?

You are very welcome to update the code in the example.