Calculate Image Moments#

Synopsis#

Calculate image moments from binary images and an interactive ellipse.

Results#

Input Image

Input image
Output
 ImageMomentsCalculator (000002037138BE90)
   RTTI typeinfo:   class itk::ImageMomentsCalculator<class itk::Image<unsigned short,2> >
   Reference Count: 1
   Modified Time: 249
   Debug: Off
   Object Name:
   Observers:
     none
   Image: 000002037138FDD0
   Valid: 1
   Zeroth Moment about origin: 153
   First Moment about origin: [30, 50]
   Second Moment about origin: 15.0458 -8.8366
 -8.8366 15.0458

   Center of Gravity: [30, 50]
   Second central moments: 15.0458 -8.8366
 -8.8366 15.0458

   Principal Moments: [950, 3654]
   Principal axes: 0.707107 0.707107
 -0.707107 0.707107

Jupyter Notebook#

https://mybinder.org/badge_logo.svg

Code#

Python#

#!/usr/bin/env python

import sys
import itk
import argparse

parser = argparse.ArgumentParser(description="Calculate Image Moments.")
parser.add_argument("input_image")

args = parser.parse_args()
input_image = itk.imread(args.input_image)
momentor = itk.ImageMomentsCalculator.New(Image=input_image)
momentor.Compute()
print(momentor)

C++#

#include "itkImageFileReader.h"
#include "itkImageMomentsCalculator.h"

int
main(int argc, char * argv[])
{
  if (argc != 2)
  {
    std::cerr << "Usage: " << std::endl;
    std::cerr << argv[0];
    std::cerr << " <InputFileName>" << std::endl;
    std::cerr << "Prints image moments from unsigned short binary image.";
    std::cerr << std::endl;
    return EXIT_FAILURE;
  }
  const char * inputFileName = argv[1];

  constexpr unsigned int Dimension = 2;

  using PixelType = unsigned short;
  using ImageType = itk::Image<PixelType, Dimension>;

  const auto input = itk::ReadImage<ImageType>(inputFileName);

  // clang-format off
  using FilterType = itk::ImageMomentsCalculator<ImageType>;
  // clang-format on
  auto filter = FilterType::New();
  filter->SetImage(input);
  filter->Compute();
  filter->Print(std::cout, 0);

  return EXIT_SUCCESS;
}

Classes demonstrated#

template<typename TImage>
class ImageMomentsCalculator : public itk::Object

Compute moments of an n-dimensional image.

This class provides methods for computing the moments and related properties of a single-echo image. Computing the (non-central) moments of a large image can easily take a million times longer than computing the various other values derived from them, so we compute the moments only on explicit request, and save their values (in an ImageMomentsCalculator object) for later retrieval by the user.

The non-central moments computed by this class are not really intended for general use and are therefore in index coordinates; that is, we pretend that the index that selects a particular pixel also equals its physical coordinates. The center of gravity, central moments, principal moments and principal axes are all more generally useful and are computed in the physical coordinates defined by the Origin and Spacing parameters of the image.

The methods that return values return the values themselves rather than references because the cost is small compared to the cost of computing the moments and doing so simplifies memory management for the caller.

See itk::ImageMomentsCalculator for additional documentation.