Skip to content

GenerateModelsFromLabels

Repository source: GenerateModelsFromLabels

Description

This example uses vtkDiscreteFlyingEdges3D or vtkDiscreteMarchingCubes to create vtkPolyData models from a 3D volume that contains discrete labels. These volumes are normally the output of a segmentation algorithm. The polydata for each label will be output into a separate file.

You can load these files into ParaView, where they will appear as a series of time steps. You can then single step through displaying the polydate from each file making up the series.

If you want to see the segmentation results as cube models, see the example GenerateCubesFromLabels

The input volume must be in MetaIO format.

Once you generate the models, you can view them with Paraview

Note

This original source code for this example is here.

Other languages

See (Python), (PythonicAPI)

Question

If you have a question about this example, please use the VTK Discourse Forum

Code

GenerateModelsFromLabels.cxx

//
// GenerateModelsFromLabels
//   Usage: GenerateModelsFromLabels InputVolume Startlabel Endlabel
//          where
//          InputVolume is a meta file containing a 3 volume of
//            discrete labels.
//          StartLabel is the first label to be processed
//          EndLabel is the last label to be processed
//          NOTE: There can be gaps in the labeling. If a label does
//          not exist in the volume, it will be skipped.
//
//
#include <vtkDiscreteFlyingEdges3D.h>
#include <vtkDiscreteMarchingCubes.h>
#include <vtkGeometryFilter.h>
#include <vtkImageAccumulate.h>
#include <vtkImageData.h>
#include <vtkMaskFields.h>
#include <vtkMetaImageReader.h>
#include <vtkNew.h>
#include <vtkPointData.h>
#include <vtkThreshold.h>
#include <vtkWindowedSincPolyDataFilter.h>
#include <vtkXMLPolyDataWriter.h>

#include <vtk_cli11.h>
#include <vtk_fmt.h>
// clang-format off
#include VTK_FMT(fmt/format.h)
// clang-format on

#include <iostream>

namespace fs = std::filesystem;

int main(int argc, char* argv[])
{
  CLI::App app{
      "Creates vtkPolyData models from a 3D volume that contains discrete "
      "labels. "
      "These volumes are normally the output of a segmentation algorithm."
      " The polydata for each label will be output into a separate file."};

  // Define options
  std::string fileName;
  app.add_option("fileName", fileName, "Input volume e.g. Frog/frogtissue.mhd.")
      ->required()
      ->check(CLI::ExistingFile);
  unsigned int startLabel{1};
  app.add_option("startLabel", startLabel,
                 "The starting label in the input volume, e.g. 1.");
  unsigned int endLabel{29};
  app.add_option("endLabel", endLabel,
                 "The ending label in the input volume, e.g. 29.");
  bool flyingEdges = true;
  app.add_flag("-m{false},!-n", flyingEdges,
               "Use flying edges by default, marching cubes if set.");

  CLI11_PARSE(app, argc, argv);

  if (flyingEdges)
  {
#define USE_FLYING_EDGES
  }

  if (startLabel > endLabel)
  {
    auto tmp = endLabel;
    endLabel = startLabel;
    startLabel = tmp;
  }

  // Create all of the classes we will need.
  vtkNew<vtkMetaImageReader> reader;
  vtkNew<vtkImageAccumulate> histogram;
#ifdef USE_FLYING_EDGES
  vtkNew<vtkDiscreteFlyingEdges3D> discreteCubes;
#else
  vtkNew<vtkDiscreteMarchingCubes> discreteCubes;
#endif
  vtkNew<vtkWindowedSincPolyDataFilter> smoother;
  vtkNew<vtkThreshold> selector;
  vtkNew<vtkMaskFields> scalarsOff;
  vtkNew<vtkGeometryFilter> geometry;
  vtkNew<vtkXMLPolyDataWriter> writer;

  // Define all of the variables
  std::string filePrefix = "Label";
  unsigned int smoothingIterations = 15;
  double passBand = 0.001;
  double featureAngle = 120.0;

  // Generate models from labels.
  // 1) Read the meta file
  // 2) Generate a histogram of the labels
  // 3) Generate models from the labeled volume
  // 4) Smooth the models
  // 5) Output each model into a separate file

  reader->SetFileName(fileName.c_str());

  histogram->SetInputConnection(reader->GetOutputPort());
  histogram->SetComponentExtent(0, endLabel, 0, 0, 0, 0);
  histogram->SetComponentOrigin(0, 0, 0);
  histogram->SetComponentSpacing(1, 1, 1);
  histogram->Update();

  discreteCubes->SetInputConnection(reader->GetOutputPort());
  discreteCubes->GenerateValues(endLabel - startLabel + 1, startLabel,
                                endLabel);

  smoother->SetInputConnection(discreteCubes->GetOutputPort());
  smoother->SetNumberOfIterations(smoothingIterations);
  smoother->BoundarySmoothingOff();
  smoother->FeatureEdgeSmoothingOff();
  smoother->SetFeatureAngle(featureAngle);
  smoother->SetPassBand(passBand);
  smoother->NonManifoldSmoothingOn();
  smoother->NormalizeCoordinatesOn();
  smoother->Update();

  selector->SetInputConnection(smoother->GetOutputPort());
#ifdef USE_FLYING_EDGES
  selector->SetInputArrayToProcess(0, 0, 0,
                                   vtkDataObject::FIELD_ASSOCIATION_POINTS,
                                   vtkDataSetAttributes::SCALARS);
#else
  selector->SetInputArrayToProcess(0, 0, 0,
                                   vtkDataObject::FIELD_ASSOCIATION_CELLS,
                                   vtkDataSetAttributes::SCALARS);
#endif

  // Strip the scalars from the output
  scalarsOff->SetInputConnection(selector->GetOutputPort());
  scalarsOff->CopyAttributeOff(vtkMaskFields::POINT_DATA,
                               vtkDataSetAttributes::SCALARS);
  scalarsOff->CopyAttributeOff(vtkMaskFields::CELL_DATA,
                               vtkDataSetAttributes::SCALARS);

  geometry->SetInputConnection(scalarsOff->GetOutputPort());

  writer->SetInputConnection(geometry->GetOutputPort());

  auto appFn = fs::path((app.get_name())).stem().string();
  // Write the files into the current file path.
  auto filePath = fs::current_path();
  std::cout << appFn
            << " is writing the files to: " << fs::current_path().string()
            << std::endl;
  for (unsigned int i = startLabel; i <= endLabel; i++)
  {
    // See if the label exists, if not skip it.
    double frequency =
        histogram->GetOutput()->GetPointData()->GetScalars()->GetTuple1(i);
    if (frequency == 0.0)
    {
      continue;
    }

    // Select the cells for a given label.
    selector->SetLowerThreshold(i);
    selector->SetUpperThreshold(i);

    // Output the polydata.
    auto fn = fmt::format("{:s}{:d}.vtp", filePrefix, i);
    std::cout << fn << std::endl;

    writer->SetFileName((filePath / fn).string().c_str());
    writer->Write();
  }
  return EXIT_SUCCESS;
}

CMakeLists.txt

cmake_minimum_required(VERSION 3.12 FATAL_ERROR)

project(GenerateModelsFromLabels)

find_package(VTK COMPONENTS 
  CommonCore
  CommonDataModel
  FiltersCore
  FiltersGeneral
  FiltersGeometry
  IOImage
  IOXML
  ImagingStatistics
)

if (NOT VTK_FOUND)
  message(FATAL_ERROR "GenerateModelsFromLabels: Unable to find the VTK build folder.")
endif()

# Prevent a "command line is too long" failure in Windows.
set(CMAKE_NINJA_FORCE_RESPONSE_FILE "ON" CACHE BOOL "Force Ninja to use response files.")
add_executable(GenerateModelsFromLabels MACOSX_BUNDLE GenerateModelsFromLabels.cxx )
  target_link_libraries(GenerateModelsFromLabels PRIVATE ${VTK_LIBRARIES}
)
# vtk_module_autoinit is needed
vtk_module_autoinit(
  TARGETS GenerateModelsFromLabels
  MODULES ${VTK_LIBRARIES}
)

Download and Build GenerateModelsFromLabels

Click here to download GenerateModelsFromLabels and its CMakeLists.txt file. Once the tarball GenerateModelsFromLabels.tar has been downloaded and extracted,

cd GenerateModelsFromLabels/build

If VTK is installed:

cmake ..

If VTK is not installed but compiled on your system, you will need to specify the path to your VTK build:

cmake -DVTK_DIR:PATH=/home/me/vtk_build ..

Build the project:

make

and run it:

./GenerateModelsFromLabels

WINDOWS USERS

Be sure to add the VTK bin directory to your path. This will resolve the VTK dll's at run time.