New Functionalities for Spherical Demons Registration
ITK has recently added the ability to register surfaces in the spherical domain based on the quad edge mesh data structure. This includes both rotational and deformable registrations of two meshes [1, 2]. The deformable registration of two meshes was implemented using the spherical demons registration as proposed by Yeo et al. [3]. Based on this work, it is now possible to perform automated atlas based labeling of surfaces generated by the nonlinear mapping between an atlas surface and a subject-specific surface. We have developed several ITK filters that will facilitate such a mapping and provide the ability to evaluate such atlas-based registration schemes.
Scalar Assignment
itk::AssignScalarsMeshFilter is a filter that iterates through PointDataContainer and assigns each scalar of the PointData from one mesh to another [4]. It takes two meshes, one as an input and the other as a source, to assign scalar values of the nodes from the source mesh to the input mesh. Both of the two meshes should have the same number of nodes. The filter can be used in the case where the user needs to pass scalar values from one mesh to another mesh. For example, when the mesh is deformed in the process of registration, this filter can be used to assign any scalar value from the initial mesh to the deformed mesh to visualize the effects of the deformation onto scalar values on the mesh. During non-linear registration, we have used this filter to assign anatomical labels based on a warped atlas surface onto the current subject surface. The following lines show how to use the AssignScalarsMeshFilte
| typedef itk::AssignScalarsMeshFilter<meshType, meshType> AssignFilterType; AssignFilterType::Pointer assignFilter = AssignFilterType::New(); assignFilter-> SetInputMesh(fixedMeshReader->GetOutput()); assignFilter-> SetSourceMesh(demonsFilter->GetOutput()); assignFilter->Update( ); |
Rescale Scalars
The simple but useful mesh filter, itk::RescaleScalarsMeshFiler, does linear transformation of the scalar values assigned to a mesh [5]. It is an extension of itk::RescaleIntensityImageFilter [6] from images. The minimum and maximum scalar values of the input mesh are calculated automatically and then a point-wise linear transformation is performed on the input scalar values by using the output minimum and maximum that are specified by the user as shown here.
| typedef itk:: RescaleScalarsMeshFilter <meshType, meshType> RescaleFilterType; RescaleFilterType::Pointer rescaleFilter = RescaleFilterType::New(); rescaleFilter-> SetInput(fixedMeshReader->GetOutput()); rescaleFilter->SetOutputMinimum(outputMinimum); rescaleFilter-> SetOutputMaximum(outputMaximum); rescaleFilter->Update( ); |
Mesh To List Adaptor
The itk::MeshToListAdaptor filter can generate a list of vectors from the scalars assigned to the input mesh [7]. It provides the necessary conversion allowing the user to utilize the itk::Statistics::SampleToHistogramFilter [6] for histogram generation of mesh scalar values. Here is the example on how to use this filter to generate a ListSample from the input mesh.
| typedef itk::MeshToListAdaptor<meshType, ListSampleType> MeshAdapterType; MeshAdapterType::Pointer meshtoListFilter = MeshAdapterType::New(); meshtoListFilter->SetMesh(MeshReader->GetOutput()); meshtoListFilter->Compute( ); |
In this filter, the Compute() method is used to trigger the generation of the ListSample. To generate the histogram, the itk::Statistics::SampleToHistogramFilter can be used to calculate the frequency of each bin on the histogram. Figure1(a) shows a spherical mesh with assigned scalars between 0 and 1. Figure 1(b) shows the histogram resulting from itk::Statistics::SampleToHistogramFilter.
Figure 1: (a) Spherical mesh with scalars assigned between 0.0 and 1.0. (b) Resulting histogram generated
from itk::Statistics::SampleToHistogramFilter. The resulting histogram was generated using 128 bins.
Histogram Matching
This filter, itk::HistogramMatchingMeshFilter, normalizes the scalar values of an input mesh based on the scalar values of a reference mesh [8]. It is an extension of the itk::HistogramMatchingImageFilter. This filter can be used to match the PointData from the moving mesh to the fixed mesh based on the histogram. The filter takes the source mesh and reference mesh as its inputs. The scalar values of the source mesh will be changed to match the histogram for scalar values from the reference mesh.
| typedef itk::HistogramMatchingQuadEdgeMeshFilter < MeshType, MeshType > FilterType; FilterType::Pointer filter = FilterType::New(); filter->SetSourceMesh(srcReader->GetOutput()); filter->SetReferenceMesh(refReader->GetOutput()); |
SetNumberOfHistogramLevels() sets the number of bins used when creating histograms of the source and reference meshes. SetNumberOfMatchPoints() governs the number of quantile values to be matched.
| filter->SetNumberOfHistogramLevels(1024); filter->SetNumberOfMatchPoints(7); filter->Update( ); |
Figure 2 demonstrates the effectiveness of the histogram matching mesh filter on two simple meshes as shown in Figure 2(a) and Figure 2(b), along with their scalar values in the range of [0.0, 1.0] and [0.5, 1.0] respectively. The output mesh with histogram matched scalars is shown in Figure 3(c).
Figure 2: The effect of histogram matching mesh filter. The scalar values of meshes in this figure are colored by
the sidebar shown on the right. (a) the input/source mesh of the filter. (b) the reference mesh of the filter.
(c) the output of the filter.
Warping Meshes
The itk::WarpMeshFilter warps a mesh using an input deformation field [9]. It was originally designed to warp anatomical labels that are mapped from the atlas space onto a subject surface. The WarpMeshFilter takes three inputs: 1) input mesh, 2) reference mesh and 3) deformation field. The input mesh contains the information of initial point positions while the reference mesh contains the scalar information that is going to be interpolated. The points on the reference mesh could be the same as the input mesh or not. The deformation field can be provided as a mesh with vectors (instead of scalars) associated with its points. The input mesh and deformation field should have the point-to-point correspondence with each other.
Each point on the input mesh is moved according to the corresponding vector on the deformation field (pnew= pin + d). A new scalar value is calculated by using interpolation on the reference mesh and mapped back onto the input mesh. An interpolator is needed by the filter since pnew is likely to fall inside of a cell and not directly on a vertex. Figure 3 shows the results of applying the warp mesh filter.
Figure 3: Shows (a) the input fixed mesh and the reference moving mesh. (b) the deformation field on the mesh. Vectors are shown as arrows on the mesh with the magnitudes on the associated color map. (c) the output of the filter. The scalar values on (a) and (c) are between 0.0 and 1.0 and colored from blue to red.
Mesh Similarity Calculator
The itk::MeshSimilarityCalculator filter is used to calculate the similarity of two labeled meshes [10]. The input meshes are assumed to have labels assigned to the surface as scalar integral values and the surfaces have the same geometry. The filter supports both the Dice and Jaccard indexes.
The actual overlap of A and B here are calculated based on the surface area of each label. If all of the points of a cell are associated with the same label, the whole triangular area of the cell is added to the calculation of A or B. However, if only one point of that cell is associated with the label, a third of the triangle area is added. If there are two points associated with the label, then two thirds of the area is added.
In order to calculate the Dice or Jaccard index, the value of the label that needs to be evaluated at both of the input meshes should be specified by the user.
| typedef itk::QuadEdgeMeshSimilarityCalculator < InputMeshType1, InputMeshType2 > SimilarityCalculatorType; SimilarityCalculatorType::Pointer similarityCalculator = SimilarityCalculatorType::New(); similarityCalculator->SetLabelValue(1); similarityCalculator->Update(); |
This filter can be used to evaluate the overlapping of warped labels on the subject’s cortical surface with the manual labels delineated by raters.
Conclusions and Future Work
The presented filters were originally created for the purpose of developing an automated cortical labeling algorithm using ITK. However, many of these filters could be used in several applications. These filters are being added to a new major update to ITK (version 4). This will provide users with the necessary infrastructure to develop novel, based surface registration algorithms using ITK.
Acknowledgements
This work is funded by NIH/NINDS grant NS050568 and thanks to Luis Ibáñez for his help and discussion.
Reference
[1] L. Ibanez, M. Audette, B.T.T. Yeo, and P. Golland. Rotational Registration of Spherical Surfaces Represented as
QuadEdge Meshes. Insight Journal. 2009.
[2] L. Ibanez, M. Audette, B.T.T. Yeo, and P. Golland. Spherical Demons Registration of Spherical Surfaces.
Insight Journal. 2009.
[3] B. T. T. Yeo, M. Sabuncu, T. Vercauteren, N. Ayache, B. Fisch, and P. Golland.
Spherical demons: Fast diffeomorphic landmark-free surface registration. IEEE TMI, 29(3): 650 – 668, 2009.
[4] W. Li and V.A. Magnotta. Assign Scalars Mesh Filter. Insight Journal. 2010.
[5] W. Li and V.A. Magnotta. Rescale Scalars Mesh Filter. Insight Journal. 2010.
[6] L. Ibanez, W. Schroeder, L. Ng, and J. Cates. The ITK Software Guide. Kitware, Inc. ISBN 1-930934-15-7,
http://www.itk.org/ItkSoftwareGuide.pdf, second edition, 2005.
[7] W. Li, V.A. Magnotta and L. Ibanez. Mesh To List Adaptor. Insight Journal. 2010.
[8] W. Li and V.A. Magnotta. Histogram Matching Mesh Filter. Insight Journal. 2010.
[9] W. Li and V.A. Magnotta. Warp Mesh Filter. Insight Journal. 2010.
[10] W. Li and V.A. Magnotta. Mesh Similarity Calculator. Insight Journal. 2010.
Wen Li is currently pursuing her PhD in the department of Biomedical Engineering at the University of Iowa and is employed as a graduate research assistant in the department of Radiology. She got her B.S. and M.S. degrees at Shanghai Jiao Tong University in China and spent three years as a lecturer and research assistant at the Science and Technology University of Shanghai. Her current research is automated surface parcellation of the human cerebral cortex using MR images.
Vincent Magnotta is an Associate Professor of Radiology, Biomedical Engineering, and Psychiatry at the University of Iowa. He is currently a developer on the ITKv4 refactoring. He is also working on two open-source projects that build upon ITK and VTK: BRAINS (http://www.nitrc.org/projects/brains) and IA-FEMesh (http://www.ccad.uiowa.edu/mimx/IA-FEMesh/).