[Insight-users] I am learning active contour algorithm , there is problem about itkGeodesicActiveContourLevelSetImageFilter

Bill Lorensen bill.lorensen at gmail.com
Sun Aug 11 07:52:22 EDT 2013


Be sure to build ITK and your example Release and not Debug.



On Sun, Aug 11, 2013 at 4:40 AM, zhq <15891495523 at 126.com> wrote:

> Hello , Dear:
>      I am learning itkGeodesicActiveContourLevelSetImageFilter , and I can
> segment a 2 dimension picture . And , now , I am trying to segment a 3
> dimension data . But , when I segment a 2 dimension picture , it costs me
> about one minute . But , it costs me more than 20 minutes to segment a 3
> dimension data , the data consists in 4 two-dimensional pictures  . Is my
> code wrong ?
> #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"
> #include "itkGDCMImageIO.h"
> #include "itkJPEGImageIO.h"
> #include "itkChangeInformationImageFilter.h"
> #include "itkImage.h"
> #include "itkImageSeriesReader.h"
> #include "itkNumericSeriesFileNames.h"
> #include "itkImageFileWriter.h"
> int main()
> {
>   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( -80.0 );
>   thresholder->SetUpperThreshold(     2 );
>
>   thresholder->SetOutsideValue(  0  );
>   thresholder->SetInsideValue(  255 );
> //  *************  reader   ******************************  //
>
> //  typedef  itk::ImageFileReader< InternalImageType > ReaderType;
> //  ReaderType::Pointer reader = ReaderType::New();
> //  reader->SetFileName("C:\\Users\\zhq\\Desktop\\VolumeData.vtk");
> //
>  reader->SetFileName("C:\\Users\\zhq\\Desktop\\BrainProtonDensitySlice.png");
> //  reader->SetImageIO(itk::GDCMImageIO::New());
> //
>  reader->SetFileName("C:\\Users\\zhq\\Desktop\\data\\SNAP_CR\\E403434298\\E403434298S1901I301.dcm");
> //  reader->SetImageIO(itk::JPEGImageIO::New());
> //  reader->SetFileName("C:\\Users\\zhq\\Desktop\\data\\picture8\\1.jpeg");
> //  reader->SetFileName("C:\\Users\\zhq\\Desktop\\\\1.jpg");
>
>   typedef itk::ImageSeriesReader<InternalImageType> ReaderType ;
>   ReaderType::Pointer reader = ReaderType::New();
>   typedef itk::NumericSeriesFileNames NameGeneratorType ;
>   NameGeneratorType::Pointer nameGenerator = NameGeneratorType::New();
>
>
> nameGenerator->SetSeriesFormat("C:\\Users\\zhq\\Desktop\\data\\SNAP_CR\\E403434298\\E403434298S1901I%3d.dcm");
>   nameGenerator->SetStartIndex(301);
>   nameGenerator->SetEndIndex(304);
>   nameGenerator->SetIncrementIndex(1);
>
>   reader->SetImageIO(itk::GDCMImageIO::New());
>   reader->SetFileNames(nameGenerator->GetFileNames());
>   reader->Update();
>   std::cout<<"reader update"<<std::endl;
> // *************  finish reader  **********************    //
>
>
> // *************  writer  *******************************  //
>   typedef  itk::ImageFileWriter<  OutputImageType  > WriterType;
>   WriterType::Pointer writer = WriterType::New();
>   writer->SetFileName("C:\\Users\\zhq\\Desktop\\ActiveContour.vtk");
>
>
> // *************  finish  writer   ********************** //
>
>   typedef itk::RescaleIntensityImageFilter<
>                                InternalImageType,
>                                OutputImageType >   CastFilterType;
>
>   typedef   itk::GradientMagnitudeRecursiveGaussianImageFilter<
>                                InternalImageType,
>                                InternalImageType >  GradientFilterType;
>   typedef   itk::SigmoidImageFilter<
>                                InternalImageType,
>                                InternalImageType >  SigmoidFilterType;
>
>   GradientFilterType::Pointer  gradientMagnitude =
> GradientFilterType::New();
>
>   SigmoidFilterType::Pointer sigmoid = SigmoidFilterType::New();
>
>   sigmoid->SetOutputMinimum(  0.0  );
>   sigmoid->SetOutputMaximum(  1.0  );
>
>   typedef  itk::FastMarchingImageFilter<
>                               InternalImageType,
>                               InternalImageType >
>  FastMarchingFilterType;
>
>   FastMarchingFilterType::Pointer  fastMarching =
> FastMarchingFilterType::New();
>
>   typedef itk::ChangeInformationImageFilter<InternalImageType>
> ChangeInformationFilterType ;
>   ChangeInformationFilterType::Pointer changeInformationFilter =
> ChangeInformationFilterType::New();
>   changeInformationFilter->ChangeSpacingOn();
>   changeInformationFilter->ChangeOriginOn();
>   changeInformationFilter->ChangeDirectionOn();
>
>   typedef  itk::GeodesicActiveContourLevelSetImageFilter<
> InternalImageType,
>                 InternalImageType >    GeodesicActiveContourFilterType;
>   GeodesicActiveContourFilterType::Pointer geodesicActiveContour =
>
>  GeodesicActiveContourFilterType::New();
>
>   const double propagationScaling = 0.25;
>
>   geodesicActiveContour->SetPropagationScaling( propagationScaling );
>   geodesicActiveContour->SetCurvatureScaling( 1 );
>   geodesicActiveContour->SetAdvectionScaling( 0.5 );
>
>   geodesicActiveContour->SetMaximumRMSError( 0.02 );
>   geodesicActiveContour->SetNumberOfIterations( 800 );
>
>
>   changeInformationFilter->SetInput(reader->GetOutput());
>   gradientMagnitude->SetInput( changeInformationFilter->GetOutput() );
>   sigmoid->SetInput( gradientMagnitude->GetOutput() );
>
>   geodesicActiveContour->SetInput(  fastMarching->GetOutput() );
>   geodesicActiveContour->SetFeatureImage( sigmoid->GetOutput() );
>
>   thresholder->SetInput( geodesicActiveContour->GetOutput() );
>   writer->SetInput( thresholder->GetOutput() );
>
>   const double sigma = 0.55 ;
>   gradientMagnitude->SetSigma(  sigma  );
>
>   const double alpha =  -0.1;
>   const double beta  =  100;
>
>   sigmoid->SetAlpha( alpha );
>   sigmoid->SetBeta(  beta  );
>
>   typedef FastMarchingFilterType::NodeContainer  NodeContainer;
>   typedef FastMarchingFilterType::NodeType       NodeType;
>
>   NodeContainer::Pointer seeds = NodeContainer::New();
>
>   InternalImageType::IndexType  seedPosition;
>
>   seedPosition[0] = 370;//661 ; //370 ; // 81;
>   seedPosition[1] = 443;//559 ; //443 ; // 114;
>   seedPosition[2] = 1;
>
>   const double initialDistance = 5;
>
>   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 );
>
>   fastMarching->SetOutputSize(
>
>  changeInformationFilter->GetOutput()->GetBufferedRegion().GetSize() );
>   try
>   {
>   gradientMagnitude->Update();
>   std::cout<<"gradientMagnitude update"<<std::endl;
>   }
>   catch( itk::ExceptionObject & excep )
>   {
>     std::cerr << "Exception caught !" << std::endl;
>     std::cerr << excep << std::endl;
>   }
>
>   sigmoid->Update();
>   std::cout<<"sigmoid update"<<std::endl;
>
>   fastMarching->SetOutputSize(
>            reader->GetOutput()->GetBufferedRegion().GetSize() );
>
>   fastMarching->Update();
>   std::cout<<"fastMarching update"<<std::endl;
>   try
>   {
>   geodesicActiveContour->Update();
>   std::cout<<"geodesicActiveContour update"<<std::endl;
>   }
>   catch( itk::ExceptionObject & excep )
>   {
>     std::cerr << "Exception caught !" << std::endl;
>     std::cerr << excep << std::endl;
>   }
>
>   writer->Update();
>   std::cout<<"writer update"<<std::endl;
>
>
>
>   std::cout << std::endl;
>   std::cout << "Max. no. iterations: " <<
> geodesicActiveContour->GetNumberOfIterations() << std::endl;
>   std::cout << "Max. RMS error: " <<
> geodesicActiveContour->GetMaximumRMSError() << std::endl;
>   std::cout << std::endl;
>   std::cout << "No. elpased iterations: " <<
> geodesicActiveContour->GetElapsedIterations() << std::endl;
>   std::cout << "RMS change: " << geodesicActiveContour->GetRMSChange() <<
> std::endl;
>
>   return 0;
> }
>
>
>
>
> 来自网易手机号码邮箱了解更多 <http://shouji.163.com>
>
> _____________________________________
> 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.php
>
> 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
>
>


-- 
Unpaid intern in BillsBasement at noware dot com
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20130811/e9d58882/attachment.htm>


More information about the Insight-users mailing list