[Insight-users] quesition of itkShapeDetectionLevelSetImageFilter in 3D

kingaza at gmail.com kingaza at gmail.com
Mon Apr 18 10:39:25 EDT 2005


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
> 
>


More information about the Insight-users mailing list