[Insight-users] Fw: Re: point based registration using least square

Suyang Mei suyang_mei at yahoo.com
Thu May 6 17:22:52 EDT 2010


  sorry I accidentally press the wrong button to send -

---------------------------------------------------------------------
Thanks a lot, for the help!
  I did successfully use itk::LandmarkBasedTransformInitializer to generate results.
 
My question is - does this class support affine transform? looks to me
it is working fine with rigid body transform, but the results not
correct if it is affine, for example, I input the following points -

The points are 
 
X=[36,257; 77,257; 30,370; 62,370;77,323];
 
X'=[118,257;77,257;126,370;93,370;77,323];

 R =

   -0.9999    0.0110
    0.0110    0.9999

is the rotation matrix I am expecting to get (I got from MATLAB svd), but I got 

R =

   0.980838    -0.194827
   0.194827    0.980838

Below are the codes I am using .

  Did I do something wrong here?

Thanks a lot,

Suyang

/////////////////////////////////////////////////////////////
    const double nPI = 4.0 * vcl_atan( 1.0 );  

    typedef  unsigned char  PixelType;
    const unsigned int Dimension = 2;
    
    typedef itk::Image< PixelType, Dimension >  FixedImageType;
    typedef itk::Image< PixelType, Dimension >  MovingImageType;

    FixedImageType::Pointer fixedImage   = FixedImageType::New();
    MovingImageType::Pointer movingImage = MovingImageType::New();

 
   
    // Set the transform type..
    typedef itk::Rigid2DTransform< double > TransformType;
    TransformType::Pointer transform = TransformType::New();
    typedef itk::LandmarkBasedTransformInitializer< TransformType, 
            FixedImageType, MovingImageType > TransformInitializerType;
    TransformInitializerType::Pointer initializer = TransformInitializerType::New();
    initializer->DebugOn();

    // Set fixed and moving landmarks
    TransformInitializerType::LandmarkPointContainer fixedLandmarks;
    TransformInitializerType::LandmarkPointContainer movingLandmarks;
    TransformInitializerType::LandmarkPointType point;
    TransformInitializerType::LandmarkPointType tmp;

    TransformInitializerType::LandmarkPointContainer testLandmarks;
    
 
    
    //pt1
    point[0]=36.0;
    point[1]=257.0;
    fixedLandmarks.push_back(point);
    point[0]=118.0;
    point[1]=257.0;
    movingLandmarks.push_back(point);

    //pt2
    point[0]=77.0;
    point[1]=257.0;
    fixedLandmarks.push_back(point);
    point[0]=77.0;
    point[1]=257.0;
    movingLandmarks.push_back(point);

    //pt3
    point[0]=30.0;
    point[1]=370.0;
    fixedLandmarks.push_back(point);
    point[0]=126.0;
    point[1]=370.0;
    movingLandmarks.push_back(point);

    //pt4
    point[0]=62.0;
    point[1]=370.0;
    fixedLandmarks.push_back(point);
    point[0]=93.0;
    point[1]=370.0;
    movingLandmarks.push_back(point);

    //pt5
    point[0]=77.0;
    point[1]=323.0;
    fixedLandmarks.push_back(point);
    point[0]=77.0;
    point[1]=323.0;
    movingLandmarks.push_back(point);

    // test pt
    point[0]=51.0;
    point[1]=320.0;
    testLandmarks.push_back(point);

    initializer->SetFixedLandmarks(fixedLandmarks);
    initializer->SetMovingLandmarks(movingLandmarks);

    initializer->SetTransform( transform );
    initializer->InitializeTransform();

    TransformInitializerType::PointsContainerConstIterator 
      fitr = fixedLandmarks.begin();
    TransformInitializerType::PointsContainerConstIterator 
      mitr = movingLandmarks.begin();

    TransformInitializerType::PointsContainerConstIterator 
      testtr = testLandmarks.begin();

    typedef TransformInitializerType::OutputVectorType  OutputVectorType;
    OutputVectorType error;
    OutputVectorType::RealValueType tolerance = 0.1;
    bool failed = false;
   
    std::cout << transform->TransformPoint( *testtr ) << std::endl;

     std::cout << "  The transform computed was: ";
      transform->Print(std::cout);
      std::cout << " DONE " << std::endl;
 

--- On Thu, 5/6/10, Suyang Mei <suyang_mei at yahoo.com> wrote:

From: Suyang Mei <suyang_mei at yahoo.com>
Subject: Re: [Insight-users] point based registration using least square
To: "Andinet Enquobahrie" <andinet.enqu at kitware.com>
Cc: insight-users at itk.org
Date: Thursday, May 6, 2010, 2:13 PM


  Thanks a lot, for the help!
  I did successfully use itk::LandmarkBasedTransformInitializer to generate results.
  My question is - does this class support affine transform? looks to me it is working fine with rigid body transform, but the results not correct if it is affine, for example, I input the following points -

The points are 
 
X=[36,257; 
77,257; 30,370; 62,370;77,323];
 
X'=[118,257;77,257;126,370;93,370;77,323];

 R =

   -0.9999    0.0110
    0.0110    0.9999

is the rotation matrix I am expecting to get (I got from MATLAB svd), but I got 

R =

   0.980838    -0.194827
   0.194827    0.980838



--- On Wed, 5/5/10, Andinet Enquobahrie <andinet.enqu at kitware.com> wrote:

From: Andinet Enquobahrie <andinet.enqu at kitware.com>
Subject: Re: [Insight-users] point based registration using least square
To: "Suyang Mei" <suyang_mei at yahoo.com>
Cc: insight-users at itk.org
Date: Wednesday, May 5, 2010, 4:54 AM

Hi Suyang-
There is ITK landmark-based registration class that probably does what
 you are looking for. 
itk::LandmarkBasedTransformInitializer
This registration class implements an absolute orientation solution using unit quaternions. It is a closed-form solution to the least-squares problem of computing transformation parameters between two coordinate systems.

Hope that helps-Andinet
On Tue, May 4, 2010 at 8:56 PM, Suyang Mei <suyang_mei at yahoo.com> wrote:


Hi, all -

  I have a simple question here -

I'd like to perform point based image registration by using least square approach (basically finding the rotation matrix by using singular value decomposition SVD). Inputs are two sets of corresponding points, and I just need to get the transformation matrix back. 



More like the functions in Examples\Patened\ICP*.cxx, but use closed form instead of iterative approach. What is the ITK filter to perform this?

  Thanks a lot in advance.

Suyang





      
_____________________________________

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



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





-- 
Andinet Enquobahrie, Ph.D.
Kitware, Inc. - North Carolina Office
http://www.kitware.com
(919) 969-6990 x311






      


      
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20100506/9d2fd8cc/attachment-0001.htm>


More information about the Insight-users mailing list