[Insight-users] 3D Fastmarching segmentation!
Luis Ibanez
luis.ibanez at kitware.com
Mon, 12 Apr 2004 10:34:58 -0400
Hi Xujf,
You shouldn't invoke the "SetSpeedConstant()" method
if you are at the same time providing an input speed
image.
The "SetSpeedConstant()" method is intended to be used
when you don't have an input speed image and want the
front to propagate at constant speed. This is usually
done when you want to use FastMarching as a Distance
Map calculator.
Please look at the documentation of this method in
the Doxygen page of the FastMarchingImageFilter.
http://www.itk.org/Insight/Doxygen/html/classitk_1_1FastMarchingImageFilter.html
--
It is very likely that the reason why you are not
getting any segmentation is related to the image
that is being passed as input to the FastMarching
filter. Than is, the output of the Sigmoid filter.
Please save this image in a file, or connect it to
a vtkImageViewer in order to look at its values.
The speed image should look almost like a binary image,
and should have almost zero values close to the edges
where you want the segmentation surface to stop, and
should have values close to 1.0 in the regions where
you want the front to propagate continously. Note that
you can change this value of 1.0 (as the max speed) by
using the method "SetNormalizationFactor()" in the
FastMarching filter (again, look at Doxygen doc for
a description of this method).
The key for a succesful application of FastMarching
is to provide an apropriate Speed image as input. That
is the *only* source of information that FastMarching
has from your actual original image. The result of
FastMarching will be as good (or as bad) as your speed
image is.
Please look at the SoftwareGuide in order to get
a feeling on how the Speed image should look like.
In particular look at Figure 9.13 in pdf-page 372.
http://www.itk.org/ItkSoftwareGuide.pdf
Note that the Sigmoid is not the only way of computing the
Speed image. You could also use the EdgePotentialImageFilter.
http://www.itk.org/Insight/Doxygen/html/classitk_1_1EdgePotentialImageFilter.html
No matter what mechanism you use for computing
the speed image, you must take a look at this
image before running FastMarching on it.
--
If you want to play with existing implementations of
FastMarching you might want to take a look at the
demo application:
InsightApplications/
FastMarchingLevelSet
and/or
to the FastMarching plugin in VolView.
You can download the free version of VolView from
http://www.kitware.com/products/volview.html
and you will find the source code of the plugins
under:
InsightApplications/VolviewPlugins
Regards,
Luis
-------------
xujf wrote:
> Hi,all:
> I use itkFastMarchingImageFilter to implement the segmentation. I can do the segmentation well in 2D image,
> the result image is very nice. I change my 2D segmentation code a little in order to get the 3D segmentation.
> however,I can\'t see segmeanything in ntation image
> my code is :(see affix if need code in detail)
> .....................
> .....................
> itkcurfilter=itk.itkCurvatureAnisotropicDiffusionImageFilterF3F3_New()
> itkgradfilter=itk.itkGradientMagnitudeRecursiveGaussianImageFilterF3F3_New()
> itksigfilter=itk.itkSigmoidImageFilterF3F3_New()
> itkfastmarchingfilter=itk.itkFastMarchingImageFilterF3F3_New()
> .................
> .................
> seedPosition=itk.itkIndex3()
> seedPosition.SetElement(0,100)
> seedPosition.SetElement(1,100)
> seedPosition.SetElement(2,100)
>
> node=itk.itkLevelSetNodeF3()
> node.SetValue(-10)
> node.SetIndex(seedPosition)
> seeds.Initialize()
> seeds.InsertElement(0,node)
> itkfastmarchingfilter.SetTrialPoints(seeds.GetPointer())
> itkfastmarchingfilter.SetSpeedConstant(1.0) ## line 50
> itkfastmarchingfilter.SetStoppingValue(80.0) ## line 51
>
> itkcurfilter.SetInput(reader.GetOutput())
> itkgradfilter.SetInput(itkcurfilter.GetOutput())
> itksigfilter.SetInput(itkgradfilter.GetOutput())
> itkfastmarchingfilter.SetInput(itksigfilter.GetOutput())
>
> ...........................
> ...........................
> I find if I delete the line50 and line 51, the result segmetation image is
> as same as before,that means line50 and line 51 don\'t work.I don\'t know what\'s
> wrong with my code.is there any other parameter I should set in itkFastMarchingImageFilter?
>
> Thanks in advance!
> xujf
>
>
> ------------------------------------------------------------------------
>
> import sys,os
> import paths
> import InsightToolkit as itk
> from vtkMINCReader import *
> import ConnectVTKITK as connect
>
> reader=vtkMINCReader()
> vtkexport=vtkImageExport()
> reader.SetFileName("E:/atamai/2003-01-27/atamai/examplesPmw/sbrain.mnc")
>
> cast=vtkImageCast()
> cast.SetOutputScalarTypeToFloat()
> cast.SetInput(reader.GetOutput())
>
> vtkexport.SetInput(cast.GetOutput())
>
> itkimport=itk.itkVTKImageImportF3_New()
>
> connect.ConnectVTKToITKF3(vtkexport,itkimport.GetPointer())
>
> itkcurfilter=itk.itkCurvatureAnisotropicDiffusionImageFilterF3F3_New()
> itkbinfilter=itk.itkBinaryThresholdImageFilterF3US3_New()
> itkgradfilter=itk.itkGradientMagnitudeRecursiveGaussianImageFilterF3F3_New()
> itksigfilter=itk.itkSigmoidImageFilterF3F3_New()
> itkfastmarchingfilter=itk.itkFastMarchingImageFilterF3F3_New()
>
> itkcurfilter.SetTimeStep( 0.0625);
> itkcurfilter.SetNumberOfIterations( 5 );
> itkcurfilter.SetConductanceParameter( 9.0 );
>
> itksigfilter.SetOutputMinimum(0.0)
> itksigfilter.SetOutputMaximum(65535.0)
> itksigfilter.SetAlpha(-1.0)
> itksigfilter.SetBeta(5.0)
>
> itkgradfilter.SetSigma(1.0)
> seeds=itk.itkNodeContainerF3_New()
>
> seedPosition=itk.itkIndex3()
> seedPosition.SetElement(0,100)
> seedPosition.SetElement(1,100)
> seedPosition.SetElement(2,100)
>
> node=itk.itkLevelSetNodeF3()
> node.SetValue(-10)
> node.SetIndex(seedPosition)
> seeds.Initialize()
> seeds.InsertElement(0,node)
> itkfastmarchingfilter.SetTrialPoints(seeds.GetPointer())
> itkfastmarchingfilter.SetSpeedConstant(1.0)
> itkfastmarchingfilter.SetStoppingValue(80.0)
>
> itkbinfilter.SetLowerThreshold(0.0)
> itkbinfilter.SetUpperThreshold(1.0 )
> itkbinfilter.SetOutsideValue( 0 )
> itkbinfilter.SetInsideValue(65535)
>
> itkcurfilter.SetInput(itkimport.GetOutput())
> itkgradfilter.SetInput(itkcurfilter.GetOutput())
> itksigfilter.SetInput(itkgradfilter.GetOutput())
> itkfastmarchingfilter.SetInput(itksigfilter.GetOutput())
> itkfastmarchingfilter.Modified()
> itkbinfilter.SetInput(itkfastmarchingfilter.GetOutput())
>
> vtkimport=vtkImageImport()
> itkexport=itk.itkVTKImageExportF3_New()
> itkexport.SetInput(itkfastmarchingfilter.GetOutput())
>
> connect.ConnectITKF3ToVTK(itkexport.GetPointer(),vtkimport)
> #cast1=vtkImageCast()
> #cast1.SetOutputScalarTypeToFloat()
> #cast1.SetInput(vtkimport.GetOutput())
>
> lut=vtkLookupTable()
> lut.SetTableRange(0,2000)
> lut.SetSaturationRange(0,0)
> lut.SetHueRange(0,0)
> lut.SetValueRange(0,1)
>
> data=vtkimport.GetOutput()
> vtkimport.Update()
>
> filter=vtkImageDataGeometryFilter()
> filter.SetInput(data)
> filter.SetExtent(0,180,100,100,0,180)
> #print data.GetExtent()
> mapper=vtkPolyDataMapper()
> mapper.SetInput(filter.GetOutput())
> mapper.SetLookupTable(lut)
> mapper.ScalarVisibilityOn()
> mapper.SetScalarRange(0,400)
> actor=vtkActor()
> actor.SetMapper(mapper)
>
> ren = vtkRenderer()
> renWin = vtkRenderWindow()
> renWin.AddRenderer(ren)
> iren = vtkRenderWindowInteractor()
> iren.SetRenderWindow(renWin)
> ren.AddActor(actor)
>
> iren.Initialize()
> renWin.Render()
> iren.Start()
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>