[Insight-users] exception with GeodesicActiveContours 3D segmentation

[Kishore Mosaliganti]

Hi Marco,

Try setting the output image parameters of the fast marching filter:

virtual void SetOutputOrigin (OutputPointType _arg)
virtual void SetOutputRegion (OutputRegionType _arg)
virtual void SetOutputSize (const OutputSizeType &size)
virtual void SetOutputSpacing (OutputSpacingType _arg)


Kishore

On Mon, May 3, 2010 at 9:14 AM, Marco Gheza
<marcogheza4mailinglists at gmail.com> wrote:
>
> Hi,
> I'm trying to use GeodesicActiveContoursSegmentation with a 3D volume. I want to use a PET gradient image as the feature image and the CT image as the image to be segmented.
> This is my code:
> #include "itkImage.h"
> #include "itkGeodesicActiveContourLevelSetImageFilter.h"
> #include "itkCurvatureAnisotropicDiffusionImageFilter.h"
> #include "itkGradientMagnitudeRecursiveGaussianImageFilter.h"
> #include "itkSigmoidImageFilter.h"
> #include "itkFastMarchingImageFilter.h"
> #include "itkRescaleIntensityImageFilter.h"
> #include "itkBinaryThresholdImageFilter.h"
> #include "itkImageFileReader.h"
> #include "itkImageFileWriter.h"
>
> int main( int argc, char *argv[] )
> {
>   if( argc < 7 )
>     {
>     std::cerr << "Missing Parameters " << std::endl;
>     std::cerr << "Usage: " << argv[0];
>     std::cerr << " inputImageCT  inputImagePT outputImage";
>     std::cerr << " seedX seedY seedZ";
>     return 1;
>     }
>
>   typedef   float           InternalPixelType;
>   const     unsigned int    Dimension = 3;
>   typedef itk::Image< InternalPixelType, Dimension >  InternalImageType;
>   typedef unsigned char                            OutputPixelType;
>   typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
>   typedef itk::BinaryThresholdImageFilter<
>                         InternalImageType,
>                         OutputImageType    >       ThresholdingFilterType;
>
>   ThresholdingFilterType::Pointer thresholder = ThresholdingFilterType::New();
>
>   thresholder->SetLowerThreshold( -1000.0 );
>   thresholder->SetUpperThreshold(     0.0 );
>   thresholder->SetOutsideValue(  0  );
>   thresholder->SetInsideValue(  255 );
>
>   typedef  itk::ImageFileReader< InternalImageType > ReaderType;
>   typedef  itk::ImageFileWriter<  OutputImageType  > WriterType;
>   ReaderType::Pointer reader = ReaderType::New();
>   ReaderType::Pointer readerPT = ReaderType::New();
>   WriterType::Pointer writer = WriterType::New();
>   reader->SetFileName( argv[1] );
>   readerPT->SetFileName( argv[2] );
>   writer->SetFileName( argv[3] );
>
>   typedef itk::RescaleIntensityImageFilter<
>                                InternalImageType,
>                                OutputImageType >   CastFilterType;
> ///filter in order to obtain the feature image///
>   typedef   itk::CurvatureAnisotropicDiffusionImageFilter<
>                                InternalImageType,
>                                InternalImageType >  SmoothingFilterType;
>   typedef   itk::GradientMagnitudeRecursiveGaussianImageFilter<
>                                InternalImageType,
>                                InternalImageType >  GradientFilterType;
>   typedef   itk::SigmoidImageFilter<
>                                InternalImageType,
>                                InternalImageType >  SigmoidFilterType;
>   SmoothingFilterType::Pointer smoothing = SmoothingFilterType::New();
>   GradientFilterType::Pointer  gradientMagnitude = GradientFilterType::New();
>   SigmoidFilterType::Pointer sigmoid = SigmoidFilterType::New();
>   smoothing->SetTimeStep( 0.0625 );
>   smoothing->SetNumberOfIterations(  5 );
>   smoothing->SetConductanceParameter( 9.0 );
>   const double sigma = 1.0;
>   gradientMagnitude->SetSigma(  sigma  );
>   const double alpha =  -0.5;
>   const double beta  =  3.0;
>   sigmoid->SetAlpha( alpha );
>   sigmoid->SetBeta(  beta  );
>   sigmoid->SetOutputMinimum(  0.0  );
>   sigmoid->SetOutputMaximum(  1.0  );
> ///filter in order to obtain the feature image///
>
> ///level set///
>   typedef  itk::FastMarchingImageFilter<
>                               InternalImageType,
>                               InternalImageType >    FastMarchingFilterType;
>   typedef  itk::GeodesicActiveContourLevelSetImageFilter< InternalImageType,
>                 InternalImageType >    GeodesicActiveContourFilterType;
>   FastMarchingFilterType::Pointer  fastMarching = FastMarchingFilterType::New();
>   GeodesicActiveContourFilterType::Pointer geodesicActiveContour =
>                                      GeodesicActiveContourFilterType::New();
>   const double propagationScaling = 2.0;
>   geodesicActiveContour->SetPropagationScaling( propagationScaling );
>   geodesicActiveContour->SetCurvatureScaling( 1.0 );
>   geodesicActiveContour->SetAdvectionScaling( 1.0 );
>   geodesicActiveContour->SetMaximumRMSError( 0.02 );
>   geodesicActiveContour->SetNumberOfIterations( 800 );
> ///level set///
>
> ///utility for the fast marching///
>   typedef FastMarchingFilterType::NodeContainer  NodeContainer;
>   typedef FastMarchingFilterType::NodeType       NodeType;
>   NodeContainer::Pointer seeds = NodeContainer::New();
>   InternalImageType::IndexType  seedPosition;
>
>   seedPosition[0] = atoi( argv[4] );
>   seedPosition[1] = atoi( argv[5] );
>   seedPosition[1] = atoi( argv[6] );
>   const double initialDistance = 5.0;
>   NodeType node;
>   const double seedValue = - initialDistance;
>
>   node.SetValue( seedValue );
>   node.SetIndex( seedPosition );
>   seeds->Initialize();
>   seeds->InsertElement( 0, node );
>   fastMarching->SetTrialPoints(  seeds  );
>   fastMarching->SetSpeedConstant( 1.0 );
> ///utility for the fast marching///
>
>   // 1° step
>   smoothing->SetInput( readerPT->GetOutput() );
>   gradientMagnitude->SetInput( smoothing->GetOutput() );
>   sigmoid->SetInput( gradientMagnitude->GetOutput() );
>   // 2° step
>   geodesicActiveContour->SetInput(  fastMarching->GetOutput() );
>   geodesicActiveContour->SetFeatureImage( sigmoid->GetOutput() );
>   // 3° step
>   thresholder->SetInput( geodesicActiveContour->GetOutput() );
>   writer->SetInput( thresholder->GetOutput() );
>
>   try
>     {
>     writer->Update();
>     }
>   catch( itk::ExceptionObject & excep )
>     {
>     std::cerr << "Exception caught !" << std::endl;
>     std::cerr << excep << std::endl;
>     }
> }
>
> unfortunatly, when i try to run it, i obtain this output message:
> Exception caught !
> itk::ExceptionObject (009EEB20)
> Location: "__thiscall itk::ImageConstIterator<class itk::Image<class itk::FixedArray<float,3>,3> >:
> ImageConstIterator(const class itk::Image<class itk::FixedArray<float,3>,3> *,const class itk::Imag
> Region<3> &)"
> File: c:\insighttoolkit-3.16.0\code\common\itkImageConstIterator.h
> Line: 182
> Description: itk::ERROR: Region ImageRegion (009EE8DC)
>   Dimension: 3
>   Index: [0, 0, 0]
>   Size: [128, 128, 287]
>  is outside of buffered region ImageRegion (024035E0)
>   Dimension: 3
>   Index: [0, 0, 0]
>   Size: [16, 16, 16]
>
> i don't know why. Please help me, thanks, bye,
> Marco
> _____________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Kitware offers ITK Training Courses, for more information visit:
> http://www.kitware.com/products/protraining.html
>
> Please keep messages on-topic and check the ITK FAQ at:
> http://www.itk.org/Wiki/ITK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.itk.org/mailman/listinfo/insight-users
>