[Insight-users] please help. filter pipeline question

[kdsfinger at gmail.com]

hi, all
I am trying to use filter pipeline to read in a series of Dicom
images, then rescale them to 1:1:1, followed by cast to unsigned char,
and finally write as a series of png files. But all my output images
are of 0 intensity. Can someone point out what's wrong with the
pipeline in my code? Thanks for comments.

#include "itkImage.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"

#include "itkResampleImageFilter.h"

#include "itkAffineTransform.h"
#include "itkNearestNeighborInterpolateImageFunction.h"

#include "itkCastImageFilter.h"

#include "itkImageSeriesReader.h"
#include "itkImageSeriesWriter.h"
#include "itkNumericSeriesFileNames.h"
#include "itkPNGImageIO.h"

#include "itkGDCMImageIO.h"
#include "itkGDCMSeriesFileNames.h"

#include <iostream>
#include <math.h>
#include <stdlib.h>
#include <string>

int main()
{
	typedef signed short    PixelType_ss;
	typedef unsigned char PixelType_uc;
  const unsigned int      Dimension = 3;

  typedef itk::Image< PixelType_ss, Dimension >      ImageType_ss;
  typedef itk::Image< PixelType_uc, Dimension > ImageType_uc;
  typedef itk::ImageSeriesReader< ImageType_ss >     ReaderType_ss;

  typedef itk::GDCMImageIO                        ImageIOType;
  typedef itk::GDCMSeriesFileNames                NamesGeneratorType;

  ImageIOType::Pointer gdcmIO = ImageIOType::New();
  NamesGeneratorType::Pointer namesGenerator = NamesGeneratorType::New();

  char temp[50];
	sprintf( temp, "../../testdataset");
  namesGenerator->SetInputDirectory( temp );

  const ReaderType_ss::FileNamesContainer & filenames =
                            namesGenerator->GetInputFileNames();

  unsigned int numberOfFilenames =  filenames.size();
  std::cout << numberOfFilenames << std::endl;
  for(unsigned int fni = 0; fni<numberOfFilenames; fni++)
    {
    std::cout << "filename # " << fni << " = ";
    std::cout << filenames[fni] << std::endl;
    }

  ReaderType_ss::Pointer reader_ss = ReaderType_ss::New();

  reader_ss->SetImageIO( gdcmIO );
  reader_ss->SetFileNames( filenames );

  try
    {
    reader_ss->Update();
    }
  catch (itk::ExceptionObject &excp)
    {
    std::cerr << "Exception thrown while writing the image" << std::endl;
    std::cerr << excp << std::endl;
    exit(1);
    }

  ImageType_ss::Pointer image = reader_ss->GetOutput();
  ImageType_ss::IndexType pixelIndex;
	pixelIndex[0] = 256; // x position
	pixelIndex[1] = 256; // y position
	pixelIndex[2] = 44; // z position
	std::cout<<image->GetPixel(pixelIndex)<<"!!!"<<std::endl;

  typedef itk::ResampleImageFilter<ImageType_ss,ImageType_ss> FilterType_ss;
  FilterType_ss::Pointer filter_ss = FilterType_ss::New();

  typedef itk::AffineTransform< double, Dimension >  TransformType;

  TransformType::Pointer transform = TransformType::New();
  filter_ss->SetTransform( transform );

  typedef itk::NearestNeighborInterpolateImageFunction<
                       ImageType_ss, double >  InterpolatorType_ss;

  InterpolatorType_ss::Pointer interpolator_ss = InterpolatorType_ss::New();
  filter_ss->SetInterpolator( interpolator_ss );

  filter_ss->SetDefaultPixelValue( 0 );

  double InputSpacing[ Dimension ];
  double OutputSpacing[ Dimension ];
  const ImageType_ss::SpacingType& sp = image->GetSpacing();
  InputSpacing[0] = sp[0];
  InputSpacing[1] = sp[1];
  InputSpacing[2] = sp[2];
  std::cout<<InputSpacing[0]<<" "<<InputSpacing[1]<<"
"<<InputSpacing[2]<<std::endl;
  double minspace = std::min(InputSpacing[0], InputSpacing[2]);
  OutputSpacing[0] = minspace;
  OutputSpacing[1] = minspace;
  OutputSpacing[2] = minspace;

  filter_ss->SetOutputSpacing( OutputSpacing );

  double origin[ Dimension ];
  origin[0] = 0.0;  // X space coordinate of origin
  origin[1] = 0.0;  // Y space coordinate of origin
  origin[2] = 0.0;

  filter_ss->SetOutputOrigin( origin );

  ImageType_ss::SizeType InputSize =
reader_ss->GetOutput()->GetLargestPossibleRegion().GetSize();
  ImageType_ss::SizeType   OutputSize;

  OutputSize[0] = InputSize[0];
  OutputSize[1] = InputSize[1];
  OutputSize[2] = (int)InputSize[2]*InputSpacing[2]/minspace;
	
  filter_ss->SetSize( OutputSize );

  filter_ss->SetInput( reader_ss->GetOutput() );
  filter_ss->Update();

  image = filter_ss->GetOutput();
  std::cout<<image->GetPixel(pixelIndex)<<"!!!"<<std::endl;
  //check to see if all pixels are convert right (Should have none 0
intensity. But they are all 0s here, which is wrong.)
  /*
  for (unsigned int x = 0; x < 255; x++)
		for (unsigned int y = 0; y < 255; y++){
			pixelIndex[0] = x;
			pixelIndex[1] = y;
			pixelIndex[2] = 1;
			if (image->GetPixel(pixelIndex) != 0)
				std::cout<<image->GetPixel(pixelIndex)<<std::endl;
		}
  */
  typedef itk::CastImageFilter<ImageType_ss, ImageType_uc > CastFilterType;
  CastFilterType::Pointer castFilter = CastFilterType::New();
  castFilter->SetInput( filter_ss->GetOutput() );

  typedef itk::Image< unsigned char, 2 > Image2DType;
	typedef itk::ImageSeriesWriter< ImageType_uc, Image2DType > WriterType_uc;
	WriterType_uc::Pointer writer_uc = WriterType_uc::New();
	writer_uc->SetInput( castFilter->GetOutput() );
	
	typedef itk::NumericSeriesFileNames OutputNameGeneratorType;
	OutputNameGeneratorType::Pointer outputNameGenerator =
OutputNameGeneratorType::New();
	
	std::string format = "../../testdataset/png/%03d.png";
	outputNameGenerator->SetSeriesFormat( format.c_str() );
	
	try
		{
		castFilter->Update();
		}
		catch( itk::ExceptionObject & excp )
		{
		std::cerr << "Exception thrown while reading the image" << std::endl;
		std::cerr << excp << std::endl;
		}

  ImageType_uc::ConstPointer inputImage = castFilter->GetOutput();
	ImageType_uc::RegionType region = inputImage->GetLargestPossibleRegion();
	ImageType_uc::IndexType start = region.GetIndex();
	ImageType_uc::SizeType size = region.GetSize();
	
	const unsigned int firstSlice = start[2];
	const unsigned int lastSlice = start[2] + size[2] - 1;
	outputNameGenerator->SetStartIndex( firstSlice );
	outputNameGenerator->SetEndIndex( lastSlice );
	outputNameGenerator->SetIncrementIndex( 1 );
	
	writer_uc->SetFileNames( outputNameGenerator->GetFileNames() );
	try
		{
		writer_uc->Update();
		}
		catch( itk::ExceptionObject & excp )
		{
		std::cerr << "Exception thrown while reading the image" << std::endl;
		std::cerr << excp << std::endl;
		}
}