Access data in cells#

Synopsis#

Just as custom data can be associated with points in the mesh, it is also possible to associate custom data with cells. The type of the data associated with the cells can be different from the data type associated with points. By default, however, these two types are the same. The following example illustrates how to access data associated with cells. The approach is analogous to the one used to access point data.

Results#

Points = 10
Cells  = 9

Cell 0 = 0
Cell 1 = 1
Cell 2 = 4
Cell 3 = 9
Cell 4 = 16
Cell 5 = 25
Cell 6 = 36
Cell 7 = 49
Cell 8 = 64

Code#

C++#

#include "itkMesh.h"
#include "itkLineCell.h"

int
main(int, char *[])
{
  using PixelType = float;
  using MeshType = itk::Mesh<PixelType, 2>;

  using CellType = MeshType::CellType;
  using LineType = itk::LineCell<CellType>;

  auto mesh = MeshType::New();

  using PointType = MeshType::PointType;
  PointType point;

  const unsigned int numberOfPoints = 10;
  for (unsigned int id = 0; id < numberOfPoints; id++)
  {
    point[0] = static_cast<PointType::ValueType>(id);              // x
    point[1] = std::log(static_cast<double>(id) + itk::Math::eps); // y
    mesh->SetPoint(id, point);
  }

  CellType::CellAutoPointer line;
  const unsigned int        numberOfCells = numberOfPoints - 1;
  for (unsigned int cellId = 0; cellId < numberOfCells; cellId++)
  {
    line.TakeOwnership(new LineType);
    line->SetPointId(0, cellId);     // first point
    line->SetPointId(1, cellId + 1); // second point
    mesh->SetCell(cellId, line);     // insert the cell
  }

  std::cout << "Points = " << mesh->GetNumberOfPoints() << std::endl;
  std::cout << "Cells  = " << mesh->GetNumberOfCells() << std::endl << std::endl;

  // assign data to cells
  for (unsigned int cellId = 0; cellId < numberOfCells; cellId++)
  {
    mesh->SetCellData(cellId, static_cast<PixelType>(cellId * cellId));
  }

  // retrieve data from cells
  for (unsigned int cellId = 0; cellId < numberOfCells; ++cellId)
  {
    PixelType value = static_cast<PixelType>(0.0);
    mesh->GetCellData(cellId, &value);
    std::cout << "Cell " << cellId << " = " << value << std::endl;
  }

  // access via an iterator
  using CellDataIterator = MeshType::CellDataContainer::ConstIterator;
  CellDataIterator cellDataIterator = mesh->GetCellData()->Begin();
  CellDataIterator end = mesh->GetCellData()->End();

  while (cellDataIterator != end)
  {
    PixelType cellValue = cellDataIterator.Value();
    // std::cout << cellValue << std::endl; //same values as before
    ++cellDataIterator;
  }

  return EXIT_SUCCESS;
}

Python#

#!/usr/bin/env python

import sys
import itk
import numpy as np

Dimension = 2
PixelType = itk.F  # float or double

MeshType = itk.Mesh[PixelType, Dimension]
PointType = itk.Point[itk.D, Dimension]

mesh = MeshType.New()

number_of_points = 10
for point_id in range(number_of_points):
    point = [float(point_id), np.log(float(point_id) + np.finfo(float).eps)]
    mesh.SetPoint(point_id, point)

number_of_cells = number_of_points - 1

# All cells are of type LINE so creating numpy of shape [number_of_cells x 2]
cells_array = np.zeros([number_of_cells, 2], dtype=np.uint64)

for cell_id in range(number_of_cells):
    cells_array[cell_id][0] = cell_id
    cells_array[cell_id][1] = cell_id + 1

cells_vector = itk.vector_container_from_array(cells_array.flatten())

# When all cells are same use the second arguement to pass the cell type
mesh.SetCellsArray(cells_vector, itk.CommonEnums.CellGeometry_LINE_CELL)

print("Points = ", mesh.GetNumberOfPoints())
print("Cells = ", mesh.GetNumberOfCells())

# Assign data to cells

# This can also be done by setting large array in one function call
# which would be more efficient than the following approach
for cell_id in range(number_of_cells):
    mesh.SetCellData(cell_id, cell_id * cell_id)


# Retrieve data from cells
cell_data = mesh.GetCellData()

# Obtain numpy array from the vector_container
cell_data_np = itk.array_from_vector_container(cell_data)

for cell_id in range(number_of_cells):
    # Demonstrating two ways of getting the element
    # First using GetElement and second using the numpy array
    if cell_id % 2 == 0:
        print("Cell ", cell_id, " = ", cell_data.GetElement(cell_id))
    else:
        temp = cell_data_np[cell_id]
        print("Cell ", cell_id, " = ", cell_data_np[cell_id])

Classes demonstrated#

template<typename TPixelType, unsigned int VDimension = 3, typename TMeshTraits = DefaultStaticMeshTraits<TPixelType, VDimension, VDimension>> class Mesh : public itk::PointSet<TPixelType, VDimension, TMeshTraits>

Implements the N-dimensional mesh structure.

Mesh is an adaptive, evolving structure. Typically points and cells are created, with the cells referring to their defining points. If additional topological information is required, then BuildCellLinks() is called and links from the points back to the cells that use them are created. This allows implicit topological information about the faces and edges of the cells to be determined. (For example, a “face” neighbor to a cell can be determined by intersection the sets of cells that use the points defining the face. This is an inherent assumption on the manifold relationship of the cells in the mesh.) In some cases, either because the mesh is non-manifold, because we wish to explicitly store information with the faces and edges of the mesh, or because performance requirements demand that boundaries are explicitly represented (the set intersection does not need to be performed); then Mesh can be further extended by adding explicit boundary assignments.

Overview: Mesh implements the N-dimensional mesh structure for ITK. It provides an API to perform operations on points, cells, boundaries, etc., but does not tie down the underlying implementation and storage. A “MeshTraits” structure is used to define the container and identifier types that will be used to access the mesh. See DefaultStaticMeshTraits for the set of type definitions needed. All types that are defined in the “MeshTraits” structure will have duplicate type alias in the resulting mesh itself.

One of the most important parts of using this mesh is how to create cells to insert into it. The cells for the mesh take two template parameters. The first is the pixel type, and should correspond exactly to that type given to the mesh. The second is a “CellTraits” which holds a sub-set of the “MeshTraits” structure definitions, and is also a member of them. Any cell which is to be inserted to a mesh should have MeshTraits::CellTraits as its second template parameter.

Usage: Mesh has three template parameters. The first is the pixel type, or the type of data stored (optionally) with points, cells, and/or boundaries. The second is the geometric dimension of the points defining the mesh. This also limits the maximum topological dimension of the cells that can be inserted. The third template parameter is the “MeshTraits” structure controlling type information for the mesh. Most users will be happy with the defaults, and will not have to worry about this third argument.

Template parameters for Mesh:

TPixelType = The type stored as data for an entity (cell, point, or boundary).

TMeshTraits = Type information structure for the mesh.

References

No reference information is available.

ITK Sphinx Examples: