[Insight-users] Transform point
Luis Ibanez
luis.ibanez at kitware.com
Tue May 13 15:24:23 EDT 2008
Hi Cristina,
Thanks for letting us know that you found the solution to the problem.
Please let us know if you run into any other difficulty,
Regards,
Luis
--------------------------
Cristina Ciavarro wrote:
> Hi Luis,
> thank you for your helpful information. The error was in item three:
> infact, I gave the wrong input points coordinates for the ITK system.
> Thank you!!!!
> Best regards,
> Cristina
>
> > Date: Tue, 6 May 2008 09:24:18 -0400
> > From: luis.ibanez at kitware.com
> > To: cristinaciavarro at hotmail.com
> > CC: insight-users at itk.org
> > Subject: Re: [Insight-users] Transform point
> >
> >
> > Hi Cristina,
> >
> >
> > Thanks for the detailed and clear description of the problem.
> >
> > The main suspect that comes to mind is that the collection of
> > input points may not be falling inside the BSpline grid.
> >
> > When that happens, a null-deformation is applied to the points,
> > and only the Bulk transform takes effect in changing their
> > coordinates.
> >
> >
> >
> > Potential reasons for this problem may be:
> >
> >
> > 1) Lack of setting the Direction of the BSpline Grid
> > based on the direction of your fixed image
> >
> > See for example, line 290 of the CVS version of
> > Insight/Examples/Registration/DeformableRegistration12.cxx
> >
> > 290: transform->SetGridDirection( gridDirection );
> >
> >
> > 2) Accidental errors on the settings of origin and
> > spacing of the BSpline Grid.
> >
> >
> > 3) It may be that the collection of input points are
> > in a coordinate system different from the Fixed
> > image coordinate system.
> >
> > How did you collected these points ?
> >
> >
> >
> >
> > Item (1) is the most likely suspect.
> >
> > What you may want to do, is to make sure that you
> > use ITK 3.6 (or the CVS version) and you may want
> > to take a look at the BSplineDeformableTransformInitializer
> > paper in the Insight Journal:
> >
> > http://insight-journal.org/midas/handle.php?handle=1926/1338
> >
> >
> > If you use this class for initializing your BSplineDeformable
> > transform, then you will rule out items (1) and (2) from the
> > list above.
> >
> >
> > Please give it a try, and let us know what you find,
> >
> >
> > Thanks
> >
> >
> > Luis
> >
> >
> > --------------------------
> > Cristina Ciavarro wrote:
> > > Dear itk users,
> > > I have developped a code in which I give as INPUT InputPoints, a
> file that contains deformable coefficients estimated by deformable
> registration, a file that contains VersorX, versorY, versorZ,
> translationX, translationY, translationZ computed by rigid registration
> and as OUTPUT OutputPoints, thas is the new coordinates of InputPoints.
> > >
> > > The steps are:
> > > 1.Registration of an image: I write in two different text files
> rototranslation values and deformable coefficients.
> > > 2.Application of the transformation calculated from precedent
> registration to an image. The code is correct because the output image
> is identical to image of point (1).
> > > 3.Application of the transformation calculated from registration to
> a PointSet.
> > >
> > > I have a problem in (3): outputPoints are only rototranslated.
> There isn't the applcation of deformable coefficients. Infact, if I
> comment the row of code transformSpline->SetBulkTransform( bulkTransform
> ); outputPoints are identical to inputPoints.
> > >
> > > What is the problem?
> > > Thank you!
> > > Best regards,
> > > Cristina
> > >
> > > #if defined(_MSC_VER)
> > > #pragma warning ( disable : 4786 )
> > > #endif
> > > #include
> > > #include
> > > #include
> > > #include "itkPointSet.h"
> > > #include "itkTransform.h"
> > > #include "itkDefaultStaticMeshTraits.h"
> > > #include "itkMultiResolutionImageRegistrationMethod.h"
> > > #include "itkMattesMutualInformationImageToImageMetric.h"
> > > #include "itkBSplineInterpolateImageFunction.h"
> > > #include "itkImage.h"
> > > #include "itkBSplineDeformableTransform.h"
> > > #include "itkVersorRigid3DTransform.h"
> > > #include "itkCenteredTransformInitializer.h"
> > >
> > > int main( int argc, char *argv[] )
> > > {
> > > if( argc < 6 )
> > > {
> > > std::cerr << "Missing Parameters " << std::endl;
> > > std::cerr << "Usage: " << argv[0];
> > > std::cerr << " fixedPointsFile outputPointsfile ";
> > > std::cerr << " DeformCoeff RigidValue InputImage TestImage";
> > >
> > > return EXIT_FAILURE;
> > > }
> > > const unsigned int Dimension=3;
> > >
> > > typedef itk::PointSet<
> double,3,itk::DefaultStaticMeshTraits>PointSetType;
> > >
> > > PointSetType::Pointer fixedPointSet = PointSetType::New();
> > > PointSetType::Pointer outputPointSet = PointSetType::New();
> > >
> > > typedef PointSetType::PointType PointType;
> > > typedef PointSetType::PointsContainer PointsContainer;
> > >
> > > PointsContainer::Pointer fixedPointContainer = PointsContainer::New();
> > > PointsContainer::Pointer outputPointContainer = PointsContainer::New();
> > >
> > > PointType fixedPoint;
> > > PointType outputPoint;
> > >
> > > // Read the file containing coordinates of fixed points.
> > > std::ifstream fixedFile;
> > > fixedFile.open( argv[1] );
> > > if( fixedFile.fail() )
> > > {
> > > std::cerr << "Error opening points file with name : " << std::endl;
> > > std::cerr << argv[1] << std::endl;
> > > return 2;
> > > }
> > > unsigned int pointId = 0;
> > > fixedFile>> fixedPoint;
> > > while( !fixedFile.eof() )
> > > {
> > > fixedPointContainer->InsertElement( pointId, fixedPoint );
> > > fixedFile>> fixedPoint;
> > > pointId++;
> > > }
> > > fixedPointSet->SetPoints( fixedPointContainer );
> > > std::cout <<
> > > "Number of fixed Points = " << fixedPointSet->GetNumberOfPoints()
> > > << std::endl;
> > >
> > > const unsigned int ImageDimension = 3;
> > > typedef signed short PixelType;
> > >
> > > typedef itk::Image< PixelType, ImageDimension> FixedImageType;
> > > typedef itk::Image< PixelType, ImageDimension> MovingImageMarkerType;
> > >
> > > const unsigned int SpaceDimension = ImageDimension;
> > > const unsigned int SplineOrder = 3;
> > > typedef double CoordinateRepType;
> > >
> > > typedef itk::BSplineDeformableTransform<
> > > CoordinateRepType,
> > > SpaceDimension,
> > > SplineOrder> TransformSplineType;
> > >
> > >
> > > typedef itk::VersorRigid3DTransform BulkTransformType;
> > >
> > > typedef itk::CenteredTransformInitializer< BulkTransformType,
> > > FixedImageType, MovingImageMarkerType> TransformInitializerType;
> > >
> > >
> > > TransformSplineType::Pointer transformSpline =
> TransformSplineType::New();
> > > BulkTransformType::Pointer bulkTransform = BulkTransformType::New();
> > > TransformInitializerType::Pointer initializer =
> TransformInitializerType::New();
> > >
> > > typedef itk::ImageFileReader< FixedImageType > FixedImageReaderType;
> > > typedef itk::ImageFileReader< MovingImageMarkerType>
> MovingImageReaderMarkerType;
> > >
> > > FixedImageReaderType::Pointer fixedImageReader =
> FixedImageReaderType::New();
> > > MovingImageReaderMarkerType::Pointer movingImageMarkerReader =
> MovingImageReaderMarkerType::New();
> > >
> > > fixedImageReader->SetFileName( argv[5] );
> > > movingImageMarkerReader->SetFileName( argv[6] );
> > > fixedImageReader->Update();
> > > movingImageMarkerReader->Update();
> > > //fixedImageReader->GetOutput()->Print(std::cout);
> > > //movingImageMarkerReader->GetOutput()->Print(std::cout);
> > >
> > >
> > > initializer->SetTransform( bulkTransform );
> > > initializer->SetFixedImage( fixedImageReader->GetOutput() );
> > > initializer->SetMovingImage( movingImageMarkerReader->GetOutput() );
> > > initializer->GeometryOn();
> > > initializer->InitializeTransform();
> > >
> > > //Read input transform file
> > > std::cout << std::endl;
> > > std::cout << "Reading input transformfile: " << argv[4] << std::endl;
> > > std::cout << std::endl;
> > > std::ifstream finTrans( argv[4] );
> > >
> > > //Set bulk transform parameters
> > > typedef BulkTransformType::ParametersType BulkParametersType;
> > > BulkParametersType bulkParameters(
> bulkTransform->GetNumberOfParameters() );
> > > std::cout << "GetNumberOfParameters bulk "<<
> > > bulkTransform->GetNumberOfParameters()<< std::endl;
> > >
> > > finTrans>> bulkParameters[0]>> bulkParameters[1]>> bulkParameters[2];
> > > finTrans>> bulkParameters[3]>> bulkParameters[4]>> bulkParameters[5];
> > > std::cout << bulkParameters<GetOutput();
> > > FixedImageType::RegionType fixedRegion =
> fixedImage->GetBufferedRegion();
> > >
> > >
> > > typedef TransformSplineType::RegionType RegionType;
> > > RegionType bsplineRegion;
> > > RegionType::SizeType gridSizeOnImage;
> > > RegionType::SizeType gridBorderSize;
> > > RegionType::SizeType totalGridSize;
> > >
> > > gridSizeOnImage.Fill( 17 );
> > > gridSizeOnImage[2]=6;
> > > gridBorderSize.Fill( 3 ); // Border for spline order = 3 ( 1 lower,
> 2 upper )
> > > totalGridSize = gridSizeOnImage + gridBorderSize;
> > >
> > > bsplineRegion.SetSize( totalGridSize );
> > >
> > > typedef TransformSplineType::SpacingType SpacingType;
> > > SpacingType spacing = fixedImage->GetSpacing();
> > >
> > > typedef TransformSplineType::OriginType OriginType;
> > > OriginType origin = fixedImage->GetOrigin();;
> > >
> > > FixedImageType::SizeType fixedImageSize = fixedRegion.GetSize();
> > >
> > > for(unsigned int r=0; r(fixedImageSize[r] - 1) /
> > > static_cast(gridSizeOnImage[r] - 1) );
> > > origin[r] -= spacing[r];
> > > }
> > >
> > > transformSpline->SetGridSpacing( spacing );
> > > transformSpline->SetGridOrigin( origin );
> > > transformSpline->SetGridRegion( bsplineRegion );
> > >
> > > std::cout << "spacing[0]= "<< spacing[0] << std::endl;
> > > std::cout << "spacing[1]= "<< spacing[1] << std::endl;
> > > std::cout << "spacing[2]= "<< spacing[2] << std::endl;
> > >
> > > std::cout << "origin[0]= "<< origin[0] << std::endl;
> > > std::cout << "origin[1]= "<< origin[1] << std::endl;
> > > std::cout << "origin[2]= "<< origin[2] << std::endl;
> > >
> > >
> > > //Read input transform file
> > > std::cout << std::endl;
> > > std::cout << "Reading input transformfile: " << argv[3] << std::endl;
> > > std::cout << std::endl;
> > > std::ifstream trasformazioneElastica( argv[3] );
> > >
> > >
> > > //Set deformable transform parameters
> > > typedef TransformSplineType::ParametersType ParametersType;
> > > ParametersType finalParameters(
> transformSpline->GetNumberOfParameters() );
> > > std::cout << "GetNumberOfParameters "<<
> transformSpline->GetNumberOfParameters()<< std::endl;
> > >
> > > for (int i=0; i< 10800; i++)
> > > {
> > >
> > > trasformazioneElastica>> finalParameters[i];
> > > }
> > >
> > > typedef BulkTransformType::CenterType CenterType;
> > > CenterType center;
> > > center[0]=0;
> > > center[1]=0;
> > > center[2]=-235;
> > > bulkTransform->SetCenter( center );
> > > bulkTransform->SetParameters( bulkParameters );
> > >
> > > transformSpline->SetBulkTransform( bulkTransform );
> > >
> > > transformSpline->SetParameters( finalParameters );
> > >
> > > for(int i = 0; iGetNumberOfPoints(); i++)
> > > {
> > > fixedPointSet->GetPoint( i, &fixedPoint );
> > > outputPoint = transformSpline->TransformPoint( fixedPoint );
> > > outputPointContainer->InsertElement( i, outputPoint );
> > > }
> > > outputPointSet->SetPoints(outputPointContainer);
> > >
> > > //Write the new pointset into file
> > > std::ofstream outputFile;
> > > outputFile.open( argv[2]);
> > > if( outputFile.fail() )
> > > {
> > > return -1;
> > > }
> > > for(int i = 0; iGetNumberOfPoints(); i++)
> > > {
> > > outputPointSet->GetPoint( i, & outputPoint );
> > > for(int j = 0; j
> > > _________________________________________________________________
> > > Una cena tra amici? Cerca il ristorante con Live Search Maps
> > > http://maps.live.it
> > > _______________________________________________
> > > Insight-users mailing list
> > > Insight-users at itk.org
> > > http://www.itk.org/mailman/listinfo/insight-users
> > >
>
>
> ------------------------------------------------------------------------
> Live Search Maps Vorresti un pub diverso dal solito? Cercalo con Live
> Search Maps! <http://maps.live.it >
More information about the Insight-users
mailing list