[Insight-users] quesition of itkShapeDetectionLevelSetImageFilter in 3D

Luis Ibanez luis.ibanez at kitware.com
Wed Apr 27 12:07:33 EDT 2005



Hi Kingaza,

Please read the chapter on LevelSet methods from the
ITK Software Guide


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

it explains in detail how to compute speed images.



If you look at the code of the most recent CVS version of the LevelSet
examples in InsightApplications, all of them are using the Speed image
as input to the FastMarching preprocessing filter.


   Regards,


       Luis



--------------------------
kingaza at gmail.com wrote:

> Hi Luis,
> 
> Thanks for your answer.
> 
> but i have still question
> 1) I don't know how to generate such a speed image you mentioned. I
> think the image from fastMarching filter is just the speed image, and
> what i do is only to set the seeds and the constant speed.
> 
> 2)I know the Origin and the Spcing is very important, and in my first
> letter I have mentioned that the codes for origin and spacing are
> added by myself.
>              fastMarching->SetOutputOrigin( reader->GetOutput()->GetOrigin() );
>              fastMarching->SetOutputSpacing(
> reader->GetOutput()->GetSpacing() );
> but even if i did so, i could not get the correct result. 
> 
> how should i do?
> now i am reading ShapeDetectionLevelSet in InsightApplications, and it
> works some acceptable. but i am so unfamiliar with fltk... :(
> 
> any help is appreciated.
> 
> 
> Regards,
> Kingaza
> 
> On 4/24/05, Luis Ibanez <luis.ibanez at kitware.com> wrote:
> 
>>
>>Hi Kingaza,
>>
>>Thanks for posting your code.
>>
>>1) The FastMarching will give better results if you
>>   feed it with the speed image instead of using
>>   the constant speed:
>>
>>       fastMarching->SetSpeedConstant( 0.9 );
>>
>>   Examples on how to connect the speed image as input
>>   to the FastMarching filter are available in
>>
>>        Insight/Examples/Segmentation
>>
>>2) The Origin and Spacing parameters in FastMarcing
>>   are *VERY* important for this filter. You have
>>   them commented out in your code.
>>
>>Regards,
>>
>>
>>    Luis
>>
>>------------------------
>>kingaza at gmail.com wrote:
>>
>>
>>>hi Luis,
>>>
>>>thx for your kindly answer
>>>
>>>this function is invoked in my mfc application.
>>>and i get:
>>>ei = shapeDetection->GetElapsedIterations() = 800;
>>>RMSChange = shapeDetection->GetRMSChange() = 0.044165259529942602;
>>>
>>>while
>>>      shapeDetection->SetMaximumRMSError( 0.02 );
>>>      shapeDetection->SetNumberOfIterations( 800 );
>>>
>>>Regards,
>>>Kingaza
>>>
>>>
>>>
>>>
>>>*.h
>>>int ShapeDetection3D(string input, string output, int seednum, int*
>>>seed, int iterations=15, double sigma=1.0, double alpha=-0.5, double
>>>beta=3.0, double curvatureScaling=0.05, double propagationScaling=1.0,
>>>double initialDistance=5.0);
>>>
>>>*.cpp
>>>int ShapeDetection3D(string input, string output, int seednum, int*
>>>seed, int iterations,  double sigma, double alpha, double beta,
>>>double curvatureScaling, double propagationScaling, double
>>>initialDistance)
>>>{
>>>
>>>      const     unsigned int    Dimension = 3;
>>>      typedef itk::Image< InternalPixelType, Dimension >  InternalImageType;
>>>
>>>      typedef itk::Image< BinaryPixelType, 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();
>>>      WriterType::Pointer writer = WriterType::New();
>>>
>>>      reader->SetFileName( input.c_str() );
>>>      writer->SetFileName( output.c_str() );
>>>
>>>      reader->Update();
>>>
>>>      typedef itk::RescaleIntensityImageFilter<InternalImageType, OutputImageType>
>>>              CastFilterType;
>>>
>>>      typedef   itk::CurvatureAnisotropicDiffusionImageFilter<
>>>              InternalImageType,
>>>              InternalImageType >  SmoothingFilterType;
>>>
>>>      SmoothingFilterType::Pointer smoothing = SmoothingFilterType::New();
>>>
>>>      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::ShapeDetectionLevelSetImageFilter< InternalImageType,
>>>              InternalImageType >    ShapeDetectionFilterType;
>>>      ShapeDetectionFilterType::Pointer shapeDetection =
>>>ShapeDetectionFilterType::New();
>>>
>>>      smoothing->SetInput( reader->GetOutput() );
>>>      gradientMagnitude->SetInput( smoothing->GetOutput() );
>>>      sigmoid->SetInput( gradientMagnitude->GetOutput() );
>>>
>>>      shapeDetection->SetInput( fastMarching->GetOutput() );
>>>      shapeDetection->SetFeatureImage( sigmoid->GetOutput() );
>>>
>>>      thresholder->SetInput( shapeDetection->GetOutput() );
>>>
>>>      writer->SetInput( thresholder->GetOutput() );
>>>
>>>      smoothing->SetTimeStep( 0.0625 );//
>>>      smoothing->SetNumberOfIterations(  15 );//5
>>>      smoothing->SetConductanceParameter( 9.0 );
>>>
>>>      gradientMagnitude->SetSigma(  sigma  );
>>>      sigmoid->SetAlpha( alpha );
>>>      sigmoid->SetBeta(  beta  );
>>>
>>>      typedef FastMarchingFilterType::NodeContainer           NodeContainer;
>>>      typedef FastMarchingFilterType::NodeType                NodeType;
>>>      NodeContainer::Pointer seeds = NodeContainer::New();
>>>
>>>      InternalImageType::IndexType  seedPosition;
>>>
>>>      NodeType node;
>>>      const double seedValue = - initialDistance;
>>>
>>>      node.SetValue( seedValue );
>>>      seeds->Initialize();
>>>      for (int n=0; n< seednum; n++)
>>>      {
>>>              seedPosition[0] = seed[3*n];
>>>              seedPosition[1] = seed[3*n+1];
>>>              seedPosition[2] = seed[3*n+2];
>>>              node.SetIndex( seedPosition );
>>>              seeds->InsertElement( n, node );
>>>      }
>>>
>>>      fastMarching->SetTrialPoints(  seeds  );
>>>      fastMarching->SetSpeedConstant( 0.9 );//1.0
>>>
>>>      CastFilterType::Pointer caster1 = CastFilterType::New();
>>>      CastFilterType::Pointer caster2 = CastFilterType::New();
>>>      CastFilterType::Pointer caster3 = CastFilterType::New();
>>>      CastFilterType::Pointer caster4 = CastFilterType::New();
>>>
>>>      WriterType::Pointer writer1 = WriterType::New();
>>>      WriterType::Pointer writer2 = WriterType::New();
>>>      WriterType::Pointer writer3 = WriterType::New();
>>>      WriterType::Pointer writer4 = WriterType::New();
>>>
>>>      caster1->SetInput( smoothing->GetOutput() );
>>>      writer1->SetInput( caster1->GetOutput() );
>>>      writer1->SetFileName("ShapeDetectionLevelSetFilterOutput1.vtk");
>>>      caster1->SetOutputMinimum(   0 );
>>>      caster1->SetOutputMaximum( 255 );
>>>      writer1->Update();
>>>
>>>      caster2->SetInput( gradientMagnitude->GetOutput() );
>>>      writer2->SetInput( caster2->GetOutput() );
>>>      writer2->SetFileName("ShapeDetectionLevelSetFilterOutput2.vtk");
>>>      caster2->SetOutputMinimum(   0 );
>>>      caster2->SetOutputMaximum( 255 );
>>>      writer2->Update();
>>>
>>>      caster3->SetInput( sigmoid->GetOutput() );
>>>      writer3->SetInput( caster3->GetOutput() );
>>>      writer3->SetFileName("ShapeDetectionLevelSetFilterOutput3.vtk");
>>>      caster3->SetOutputMinimum(   0 );
>>>      caster3->SetOutputMaximum( 255 );
>>>      writer3->Update();
>>>
>>>      caster4->SetInput( fastMarching->GetOutput() );
>>>      writer4->SetInput( caster4->GetOutput() );
>>>      writer4->SetFileName("ShapeDetectionLevelSetFilterOutput4.vtk");
>>>      caster4->SetOutputMinimum(   0 );
>>>      caster4->SetOutputMaximum( 255 );
>>>
>>>
>>>      fastMarching->SetOutputSize(
>>>              reader->GetOutput()->GetBufferedRegion().GetSize() );
>>>      //      fastMarching->SetOutputOrigin(
>>>      //              reader->GetOutput()->GetOrigin() );
>>>      //      fastMarching->SetOutputSpacing(
>>>      //              reader->GetOutput()->GetSpacing() );
>>>
>>>
>>>      shapeDetection->SetPropagationScaling(  propagationScaling );
>>>      shapeDetection->SetCurvatureScaling( curvatureScaling );
>>>
>>>
>>>      shapeDetection->SetMaximumRMSError( 0.02 );
>>>      shapeDetection->SetNumberOfIterations( 800 );
>>>
>>>      try
>>>      {
>>>              writer->Update();
>>>      }
>>>      catch( itk::ExceptionObject & excep )
>>>      {
>>>              std::cerr << "Exception caught !" << std::endl;
>>>              std::cerr << excep << std::endl;
>>>              return -1;
>>>      }
>>>
>>>      int ni = shapeDetection->GetNumberOfIterations();
>>>      double RMS = shapeDetection->GetMaximumRMSError();
>>>      int ei = shapeDetection->GetElapsedIterations();
>>>      double RMSChange = shapeDetection->GetRMSChange();
>>>
>>>
>>>      std::cout << std::endl;
>>>      std::cout << "Max. no. iterations: " << ni << std::endl;
>>>      std::cout << "Max. RMS error: " << RMS << std::endl;
>>>      std::cout << std::endl;
>>>      std::cout << "No. elpased iterations: " << ei << std::endl;
>>>      std::cout << "RMS change: " << RMSChange << std::endl;
>>>
>>>              writer4->Update();
>>>
>>>
>>>
>>>
>>>      /*
>>>      typedef itk::ImageFileWriter< InternalImageType > InternalWriterType;
>>>
>>>      InternalWriterType::Pointer mapWriter = InternalWriterType::New();
>>>      mapWriter->SetInput( fastMarching->GetOutput() );
>>>      mapWriter->SetFileName("ShapeDetectionLevelSetFilterOutput4.mha");
>>>      mapWriter->Update();
>>>
>>>      InternalWriterType::Pointer speedWriter = InternalWriterType::New();
>>>      speedWriter->SetInput( sigmoid->GetOutput() );
>>>      speedWriter->SetFileName("ShapeDetectionLevelSetFilterOutput3.mha");
>>>      speedWriter->Update();
>>>
>>>      InternalWriterType::Pointer gradientWriter = InternalWriterType::New();
>>>      gradientWriter->SetInput( gradientMagnitude->GetOutput() );
>>>      gradientWriter->SetFileName("ShapeDetectionLevelSetFilterOutput2.mha");
>>>      gradientWriter->Update();
>>>      */
>>>      return 0;
>>>}
>>>
>>>On 4/18/05, Luis Ibanez <luis.ibanez at kitware.com> wrote:
>>>
>>>
>>>>Hi Kingaza,
>>>>
>>>>Thanks for posting the screen shots of your output.
>>>>
>>>>It seems that the speed image is not being connected
>>>>correctly to the FastMarching filter.
>>>>
>>>>Your manipulations of origin and spacing may easily
>>>>be part of the problem....
>>>>
>>>>You may also have an problem with the scaling of the
>>>>speed image. Are you using values from 0.0 to 1.0 ??
>>>>
>>>>Please post your code to the list.
>>>>
>>>> Thanks
>>>>
>>>>    Luis
>>>>
>>>>-----------------------------
>>>>kingaza at gmail.com wrote:
>>>>
>>>>
>>>>>hi all,
>>>>>
>>>>>I am so terribly puzzled to this problem.
>>>>>I have tried ShapeDetectionLevelSetFilter in a 3D case. Everything is
>>>>>ok, but for the last process. It's seems that the feature image out
>>>>
>>>>>from sigmoid filter doesn't work at all.
>>>>
>>>>>the attachment pictures are the results out from fastmarching, sigmoid
>>>>>and shapedetection filters, which i captured from paraview.
>>>>>
>>>>>eager for any tip!
>>>>>
>>>>>regards,
>>>>>kingaza
>>>>>
>>>>>btw, because of different spacing and origin, i add these codes
>>>>>      fastMarching->SetOutputOrigin(
>>>>>              reader->GetOutput()->GetOrigin() );
>>>>>      fastMarching->SetOutputSpacing(
>>>>>              reader->GetOutput()->GetSpacing() );
>>>>>after
>>>>>      fastMarching->SetOutputSize(
>>>>>              reader->GetOutput()->GetBufferedRegion().GetSize() );
>>>>>
>>>>>
>>>>>------------------------------------------------------------------------
>>>>>
>>>>>
>>>>>------------------------------------------------------------------------
>>>>>
>>>>>
>>>>>------------------------------------------------------------------------
>>>>>
>>>>>
>>>>>------------------------------------------------------------------------
>>>>>
>>>>>_______________________________________________
>>>>>Insight-users mailing list
>>>>>Insight-users at itk.org
>>>>>http://www.itk.org/mailman/listinfo/insight-users
>>>>
>>>>
>>>_______________________________________________
>>>Insight-users mailing list
>>>Insight-users at itk.org
>>>http://www.itk.org/mailman/listinfo/insight-users
>>>
>>>
>>
>>
> _______________________________________________
> Insight-users mailing list
> Insight-users at itk.org
> http://www.itk.org/mailman/listinfo/insight-users
> 
> 





More information about the Insight-users mailing list