Visualize a Static Sparse Malcolm 2D Level-Set Layers#

Synopsis#

Convert a binary mask into a sparse Malcolm level-set function and write its layers as an image, by sampling the function at every pixel. The mask comes from an Otsu threshold of the input image.

The written values are -1, 0, and 1.

Results#

Code#

C++#

#include "itkBinaryImageToLevelSetImageAdaptor.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkImageRegionIteratorWithIndex.h"
#include "itkMalcolmSparseLevelSetImage.h"
#include "itkOtsuMultipleThresholdsImageFilter.h"
#include "itkRescaleIntensityImageFilter.h"

int
main(int argc, char * argv[])
{
  if (argc != 3)
  {
    std::cerr << "Missing Arguments" << std::endl;
    std::cerr << argv[0] << std::endl;
    std::cerr << "<Input Image> <Output Image>" << std::endl;
    return EXIT_FAILURE;
  }

  constexpr unsigned int Dimension = 2;

  using InputPixelType = unsigned char;
  using InputImageType = itk::Image<InputPixelType, Dimension>;

  InputImageType::Pointer input = itk::ReadImage<InputImageType>(argv[1]);

  using LevelSetType = itk::MalcolmSparseLevelSetImage<Dimension>;

  // Generate a binary mask that will be used as initialization for the level
  // set.
  using OtsuFilterType = itk::OtsuMultipleThresholdsImageFilter<InputImageType, InputImageType>;
  auto otsu = OtsuFilterType::New();
  otsu->SetInput(input);
  otsu->SetNumberOfHistogramBins(256);
  otsu->SetNumberOfThresholds(1);

  using RescaleType = itk::RescaleIntensityImageFilter<InputImageType, InputImageType>;
  auto rescaler = RescaleType::New();
  rescaler->SetInput(otsu->GetOutput());
  rescaler->SetOutputMinimum(0);
  rescaler->SetOutputMaximum(1);
  rescaler->Update();

  // Convert the binary mask to a sparse level-set function
  using BinaryImageToLevelSetType = itk::BinaryImageToLevelSetImageAdaptor<InputImageType, LevelSetType>;
  auto adaptor = BinaryImageToLevelSetType::New();
  adaptor->SetInputImage(rescaler->GetOutput());
  adaptor->Initialize();

  LevelSetType::Pointer levelSet = adaptor->GetModifiableLevelSet();

  // Sample the level-set function on the image grid to show its layers
  using LayerImageType = itk::Image<LevelSetType::OutputType, Dimension>;
  auto layers = LayerImageType::New();
  layers->CopyInformation(input);
  layers->SetRegions(input->GetLargestPossibleRegion());
  layers->Allocate();

  itk::ImageRegionIteratorWithIndex<LayerImageType> it(layers, layers->GetLargestPossibleRegion());
  for (it.GoToBegin(); !it.IsAtEnd(); ++it)
  {
    it.Set(levelSet->Evaluate(it.GetIndex()));
  }

  try
  {
    itk::WriteImage(layers, argv[2]);
  }
  catch (const itk::ExceptionObject & error)
  {
    std::cerr << "Error: " << error << std::endl;
    return EXIT_FAILURE;
  }

  return EXIT_SUCCESS;
}

Classes demonstrated#

template<unsigned int VDimension> class MalcolmSparseLevelSetImage : public itk::LevelSetSparseImage<int8_t, VDimension>

Derived class for the Malcolm representation of level-set function.

This representation is a “sparse” level-set function, where values could only be { -1, 0, +1 } and organized into 1 layer { 0 }.

Template Parameters:: VDimension – Dimension of the input space

See itk::MalcolmSparseLevelSetImage for additional documentation.