itkregiongrow2dtest.cxx

DTMK软件开发包,此为开源软件,是一款很好的医学图像开发资源.

/*=========================================================================

  Program:   Insight Segmentation & Registration Toolkit
  Module:    $RCSfile: itkRegionGrow2DTest.cxx,v $
  Language:  C++
  Date:      $Date: 2008-06-10 13:02:48 $
  Version:   $Revision: 1.41 $

  Copyright (c) Insight Software Consortium. All rights reserved.
  See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.

     This software is distributed WITHOUT ANY WARRANTY; without even
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
     PURPOSE.  See the above copyright notices for more information.

=========================================================================*/
#if defined(_MSC_VER)
#pragma warning ( disable : 4786 )
#endif
// Insight classes
#include "itkImage.h"
#include "itkVector.h"
#include "vnl/vnl_matrix_fixed.h"
#include "itkImageRegionIterator.h"
#include "itkTextOutput.h"
#include "itkRGBPixel.h"

#include "itkKLMRegionGrowImageFilter.h"
#include "itkScalarImageToHistogramGenerator.h"

#define   NUMBANDS1           1
#define   NUMBANDS2           2
#define   NUMBANDS3           3
#define   NUMDIM1D            1
#define   NUMDIM2D            2
#define   NUMDIM3D            3
#define   NUMDIM4D            4
#define   NUMDIM5D            5

static unsigned int test_RegionGrowKLMExceptionHandling();
static unsigned int test_regiongrowKLM1D();
static unsigned int test_regiongrowKLM2D();
static unsigned int test_regiongrowKLM3D();
#ifndef _GLIBCXX_DEBUG
static unsigned int test_regiongrowKLM4D();
#endif

//
// This tests KLM region growing segmentation
//

int itkRegionGrow2DTest(int, char* [] )
{
  int pass;

  // Exception test the KLM algorithm
  pass = test_RegionGrowKLMExceptionHandling();
  if ( pass == EXIT_FAILURE ) return pass;

  // Test the KLM algorithm applied to 1D data
  pass = test_regiongrowKLM1D();
  if ( pass == EXIT_FAILURE ) return pass;

  // Test the KLM algorithm applied to 2D data
  pass = test_regiongrowKLM2D();
  if ( pass == EXIT_FAILURE ) return pass;

  // Test the KLM algorithm applied to 3D data
  pass = test_regiongrowKLM3D();
  if ( pass == EXIT_FAILURE ) return pass;

  // Test the KLM algorithm applied to 4D data
#ifndef _GLIBCXX_DEBUG
  pass = test_regiongrowKLM4D();
  if ( pass == EXIT_FAILURE ) return pass;
#endif
  return EXIT_SUCCESS;
}


unsigned int test_RegionGrowKLMExceptionHandling()
{
  itk::OutputWindow::SetInstance(itk::TextOutput::New().GetPointer());

  std::cout << "Testing exception handling" << std::endl;

  // Perform the exception handling testing on a 5D image

  // Generate the image data

  int sizeLen  = 3;

  typedef itk::Image<itk::Vector<double,NUMBANDS2>,NUMDIM5D> ImageType5D;
  ImageType5D::Pointer image5D  = ImageType5D::New();

  ImageType5D::SizeType imageSize5D;
  imageSize5D.Fill( sizeLen );

  ImageType5D::IndexType index5D;
  index5D.Fill(0);

  ImageType5D::RegionType region5D;

  region5D.SetSize( imageSize5D );
  region5D.SetIndex( index5D );

  image5D->SetLargestPossibleRegion( region5D );
  image5D->SetBufferedRegion( region5D );
  image5D->Allocate();
  itk::Vector<double,NUMBANDS2> pixel(0.0);
  image5D->FillBuffer(pixel);

  // Set the filter with valid inputs

  typedef itk::KLMRegionGrowImageFilter<ImageType5D,ImageType5D>
    KLMRegionGrowImageFilterType5D;

  KLMRegionGrowImageFilterType5D::Pointer
    exceptionTestingFilter5D = KLMRegionGrowImageFilterType5D::New();

  KLMRegionGrowImageFilterType5D::GridSizeType gridSize5D;
  gridSize5D.Fill(1);

  exceptionTestingFilter5D->SetInput(image5D);
  exceptionTestingFilter5D->SetGridSize(gridSize5D);
  exceptionTestingFilter5D->SetMaximumNumberOfRegions(2);
  exceptionTestingFilter5D->SetMaximumLambda(1000);

  std::cout << "Test error handling" << std::endl;

  bool passed;

#undef LOCAL_TEST_EXCEPTION_MACRO
#define LOCAL_TEST_EXCEPTION_MACRO( MSG, FILTER ) \
  passed = false; \
  try \
    { \
    std::cout << MSG << std::endl; \
    FILTER->Update(); \
    } \
  catch( itk::ExceptionObject& err ) \
    { \
    std::cout << "Caught expected error." << std::endl; \
    std::cout << err << std::endl; \
    FILTER->ResetPipeline(); \
    passed = true; \
    } \
  if ( !passed ) \
    { \
    std::cout << "Test FAILED" << std::endl; \
    return EXIT_FAILURE; \
    }


  // maximum number of regions must be greater than 1

  exceptionTestingFilter5D->SetMaximumNumberOfRegions(0);
  LOCAL_TEST_EXCEPTION_MACRO( "Maximum number of user specified region is 0",
    exceptionTestingFilter5D );

  exceptionTestingFilter5D->SetMaximumNumberOfRegions(1);
  LOCAL_TEST_EXCEPTION_MACRO( "Maximum number of user specified region is 1",
    exceptionTestingFilter5D );

  exceptionTestingFilter5D->SetMaximumNumberOfRegions(2);

  // size lengths must be divisible by the grid size along each dimension

  for( int idim = 0; idim < NUMDIM5D; idim++ )
    {
    gridSize5D[idim]++;
    exceptionTestingFilter5D->SetGridSize( gridSize5D );
    LOCAL_TEST_EXCEPTION_MACRO( "Invalid grid size, dimension " << idim+1 <<
                                " of " << NUMDIM5D, exceptionTestingFilter5D );
    gridSize5D[idim]--;
    }
  exceptionTestingFilter5D->SetGridSize( gridSize5D );

  // gridSize can't be 0

  gridSize5D[2] = 0;
  exceptionTestingFilter5D->SetGridSize( gridSize5D );
  LOCAL_TEST_EXCEPTION_MACRO( "Invalid grid size = 0",
    exceptionTestingFilter5D );
  gridSize5D[2] = 1;
  exceptionTestingFilter5D->SetGridSize( gridSize5D );

  // one region in there - what happens?

  gridSize5D.Fill( sizeLen );
  exceptionTestingFilter5D->SetGridSize( gridSize5D );
  LOCAL_TEST_EXCEPTION_MACRO( "One input region, grid too large",
    exceptionTestingFilter5D );

  std::cout << "Done testing exception handling" << std::endl;

  return EXIT_SUCCESS;
}


unsigned int test_regiongrowKLM1D()
{
  itk::OutputWindow::SetInstance(itk::TextOutput::New().GetPointer());

  std::cout << std::endl << "Begin testing one-dimension images" << std::endl;

  // Manually create an image

  typedef itk::Image<itk::Vector<unsigned char,NUMBANDS3>,NUMDIM1D> ImageType;
  typedef itk::Image<itk::Vector<double,NUMBANDS3>,NUMDIM1D> OutputImageType;

  ImageType::Pointer image  = ImageType::New();

  unsigned int numPixels = 100;
  unsigned int numPixelsHalf = 50;
  ImageType::SizeType imageSize;
  imageSize.Fill( numPixels );

  ImageType::IndexType index;
  index.Fill(0);

  ImageType::RegionType region;
  region.SetSize( imageSize );
  region.SetIndex( index );

  image->SetLargestPossibleRegion( region );
  image->SetBufferedRegion( region );
  image->Allocate();

  typedef ImageType::PixelType ImagePixelType;
  typedef itk::ImageRegionIterator< ImageType > ImageIterator;
  ImageIterator inIt( image, image->GetBufferedRegion() );
  ImageIterator inItEnd = inIt.End();

  typedef ImageType::PixelType::VectorType ImageData;
  ImageData   pixelData;
  unsigned int k = 0;
  while ( inIt != inItEnd ) {
    pixelData[0] = static_cast<unsigned char>( k );
    pixelData[1] = static_cast<unsigned char>( numPixels - k - 1 );
    if ( k < numPixelsHalf ) pixelData[2] = 47;
    else pixelData[2] = 247;
    inIt.Set( pixelData );
    ++inIt;
    ++k;
  }


  // FIRST TEST:
  // If lambda is 0, the number of final regions should equal
  // initial number of regions, the region labels should be consecutive

  // Set up the filter

  typedef itk::KLMRegionGrowImageFilter<ImageType,OutputImageType>
    KLMRegionGrowImageFilterType;

  KLMRegionGrowImageFilterType::Pointer KLMFilter =
    KLMRegionGrowImageFilterType::New();

  KLMRegionGrowImageFilterType::GridSizeType gridSize;
  gridSize.Fill(1);

  KLMFilter->SetInput( image );
  KLMFilter->SetGridSize(gridSize);

#undef LOCAL_TEST_EXCEPTION_MACRO
#define LOCAL_TEST_EXCEPTION_MACRO( FILTER ) \
  try \
    { \
    FILTER->Update(); \
    } \
  catch( itk::ExceptionObject& err ) \
    { \
    std::cout << "Caught unexpected error." << std::endl; \
    std::cout << err << std::endl; \
    return EXIT_FAILURE; \
    } \
  std::cout << std::endl << "Filter has been udpated" << std::endl


  std::cout << std::endl << "First test, lambda = 0" << std::endl;

  KLMFilter->SetMaximumLambda( 0 );
  int nregions = 2;
  KLMFilter->SetMaximumNumberOfRegions( nregions );

  // Kick off the Region grow function
  LOCAL_TEST_EXCEPTION_MACRO( KLMFilter );
  KLMFilter->Print(std::cout);

  if( numPixels != KLMFilter->GetNumberOfRegions() )
    {
    std::cout << "Test FAILED" << std::endl;
    return EXIT_FAILURE;
    }

  // Test the functions useful to test the region and border statistics
  // as the regions are merged. Primarily useful for debug operations and are
  // called several times, so prudent usage is advisable.
  // KLMFilter->PrintAlgorithmRegionStats();
  // KLMFilter->PrintAlgorithmBorderStats();

  std::cout << "Extracting and checking approximation image" << std::endl;
  std::cout << "Extracting and checking label image" << std::endl;

  // Make sure that the labelled image type is set to unsigned integer
  // as labels associated with different regions are always integers
  // This should return unique integer labels of the segmented regions.
  // The region labels should be consecutive integers beginning with 1.

  OutputImageType::Pointer outImage = KLMFilter->GetOutput();
  typedef itk::ImageRegionIterator< OutputImageType > OutputImageIterator;
  OutputImageIterator outIt( outImage, outImage->GetBufferedRegion() );

  typedef KLMRegionGrowImageFilterType::RegionLabelType LabelType;
  typedef itk::Image<LabelType, NUMDIM1D> LabelledImageType;
  LabelledImageType::Pointer labelledImage = KLMFilter->GetLabelledImage();

  typedef OutputImageType::PixelType::VectorType OutputImageData;
  ImageData         pixelIn;
  OutputImageData   pixelOut;

  typedef LabelledImageType::PixelType LabelledImagePixelType;
  typedef itk::ImageRegionIterator< LabelledImageType > LabelImageIterator;
  LabelImageIterator
    labelIt( labelledImage, labelledImage->GetBufferedRegion() );
  LabelType pixelLabel;
  LabelType m = 1;

  inIt.GoToBegin();
  while ( inIt != inItEnd )
    {
    pixelOut = outIt.Get();
    pixelIn  = inIt.Get();
    pixelLabel = labelIt.Get();

    if ( pixelOut[0] != pixelIn[0] ||
         pixelOut[1] != pixelIn[1] ||
         pixelOut[2] != pixelIn[2] ||
         pixelLabel  != m )
      {
      std::cout << "Test FAILED" << std::endl;
      if ( pixelOut[0] != pixelIn[0])
        {
        std::cout << "pixelOut[0]: " << pixelOut[0]
                  << " != "
                  << "pixelIn[0]: " << pixelIn[0]
                  << std::endl;
        }
      if ( pixelOut[1] != pixelIn[1])
        {
        std::cout << "pixelOut[1]: " << pixelOut[1]
                  << " != "
                  << "pixelIn[1]: " << pixelIn[1]
                  << std::endl;
        }
      if ( pixelOut[2] != pixelIn[2])
        {
        std::cout << "pixelOut[2]: " << pixelOut[2]
                  << " != "
                  << "pixelIn[2]: " << pixelIn[2]
                  << std::endl;
        }
      if ( pixelLabel  != m )
        {
        std::cout << "pixelLabel: " << pixelLabel
                  << " != "
                  << "m: " << m
                  << std::endl;
        }
      return EXIT_FAILURE;
      }

    ++outIt;
    ++inIt;
    ++labelIt;
    ++m;

    } //end while iterator loop

  std::cout << "Test PASSED" << std::endl;

  KLMFilter->ResetPipeline();



  // SECOND TEST:
  // merge as much as possible

  std::cout << std::endl << "Second test, merge to " << nregions <<
    " regions" << std::endl;