[Insight-users] extract filter problem

[Luis Ibanez]

Hi Michael



Please do the following:



A) Check on image size

1) In the try/catch block

 > 	try {
 > 		extractFilter->Update();
 > 	} catch (itk::ExceptionObject &e) {

replace "extractFilter" with the gradientFilter
"preFilter". Then just after Update() add

   preFilter->GetOutput()->Print( std::cout );

2) Rerun your program

3) Let us know the region sized that are printed.




B) Check for a bug in the ExtractImageFilter

   1) Replace the ExtractImageFilter with the
      "RegionOfInterestImageFilter"

   This filter is described in the SoftwareGuide

        http://www.itk.org/ItkSoftwareGuide.pdf

   Section 7.4, pdf-page 225.


If your program works with this filter, then there
is probably a buggy condition in the ExtractImageFilter.



   Please let us know what you find.


          Thanks


             Luis



---------------------------
michakuhn at gmx.ch wrote:

> Hi,
> 
> I've written the following program, which is supposed to create the gradient
> image for a particular region of the input image.
> 
> ----
> 
> #include <iostream>
> 
> #include "itkImage.h"
> #include "itkImageFileReader.h"
> #include "itkImageFileWriter.h"
> #include "itkGradientImageFilter.h"
> #include "itkGradientRecursiveGaussianImageFilter.h"
> #include "itkNumericTraits.h"
> #include "itkImageRegionConstIterator.h"
> #include "itkDiscreteGaussianImageFilter.h"
> #include "itkCastImageFilter.h"
> #include "itkExtractImageFilter.h"
> 
> int main(int argc, char** argv)
> {
> 	const unsigned int Dimension = 3;
> 	typedef short PixelType;
> 	typedef itk::Image<PixelType, Dimension> ImageType;
> 
> 	typedef double RealType;
> 
> 	typedef itk::Image<RealType, Dimension> IntermediateImageType;
> 
> 	typedef itk::CovariantVector<RealType, Dimension> GradientPixelType;
> 	typedef itk::Image<GradientPixelType, Dimension> GradientImageType;
> 	
> 	
> 	typedef itk::ImageFileReader<ImageType> ReaderType;
> 
> 	ReaderType::Pointer reader = ReaderType::New();
> 
> 	if (argc < 2) {
> 		std::cout << "image in required" << std::endl;
> 		return -1;
> 	}
> 	reader->SetFileName(argv[1]);
> 
> 	try {
> 		reader->Update();
> 	} catch (itk::ExceptionObject &e) {
> 		std::cout << e << std::endl;
> 		return -1;
> 	}
> 
> 	typedef itk::CastImageFilter<ImageType, IntermediateImageType>
> 		CastFilterType;
> 
> 	CastFilterType::Pointer castFilter = CastFilterType::New();
> 
> 	castFilter->SetInput(reader->GetOutput());
> 
> 
> 	typedef itk::DiscreteGaussianImageFilter<IntermediateImageType,
> 		IntermediateImageType>	PreFilterType;
> 
> 	PreFilterType::Pointer preFilter = PreFilterType::New();
> 	preFilter->SetInput(castFilter->GetOutput());
> 	preFilter->SetUseImageSpacing(false);
> 	preFilter->SetMaximumKernelWidth(8);
> 	std::cout << "default variance: " << preFilter->GetVariance() << std::endl;
> 	preFilter->SetVariance(1);
> 
> 	typedef itk::GradientImageFilter<IntermediateImageType, RealType, 
> 		RealType> Filter1Type;
> 
> 	Filter1Type::Pointer filter1 = Filter1Type::New();
> 	filter1->SetInput(preFilter->GetOutput());
> 
> 	typedef itk::ExtractImageFilter<GradientImageType, GradientImageType>
> 		ExtractFilterType;
> 	ExtractFilterType::Pointer extractFilter = ExtractFilterType::New();
> 
> 	extractFilter->SetInput(filter1->GetOutput());
> 
> 	typedef GradientImageType::RegionType RegionType;
> 	RegionType region;
> 	RegionType::SizeType size;
> 	RegionType::IndexType index;
> 
> 	size[0] = 1;
> 	size[1] = 32;
> 	size[2] = 32;
> 
> 	index[0] = 64;
> 	index[1] = 0;
> 	index[2] = 0;
> 
> 	region.SetSize(size);
> 	region.SetIndex(index);
> 
> 	extractFilter->SetExtractionRegion(region);
> 
> 	try {
> 		extractFilter->Update();
> 	} catch (itk::ExceptionObject &e) {
> 		std::cout << e << std::endl;
> 		return -1;
> 	}
> 
> 	GradientImageType::Pointer image1 = extractFilter->GetOutput();
> 
> 	std::cout << "image1: " << image1 << std::endl;
> 
> 	return 0;
> }
> 
> ----
> 
> This program crashes, if my input image has a size of 65x65x127 pixels. It
> works fine, however, as soon as I change (in the source code) "index[0] =
> 64" to "index[0] = 63". It still works when increasing size[0] to 2 after
> changing the index. Can somebody explain me, why I am not allowed to extract
> this particular region I've set in the code? Or is there something else
> wrong?
> 
> Thanks,
> 
> Michael
> 
> P.S.: Below the input image file (in.mhd):
> 
> ----
> 
> ObjectType = Image
> NDims = 3
> BinaryData = True
> BinaryDataByteOrderMSB = False
> ElementSpacing = 0.127664 0.127664 0.127664
> DimSize = 65 65 127
> ElementType = MET_SHORT
> ElementDataFile = in.raw
> 
> ----
>