[Insight-users] PointSet to PointSet registration problem with LevenbergMarquardtOptimizer and VersorRigid3DTransform

[Luis Ibanez]

Hi Neuhaus,

In order to make a good selection of scaling parameters
for a Versor3DTransform, you should answer the following
two questions:


A) How large is the translation that you expect,
     will be needed for registering the points ?

     1mm ?  10mm ?   2m ?


B) How large of a rotation do you anticipate to be
     necessary for registering the two datasets ?

     1 degree ? 10 degrees ?  90 degrees ?

    BTW: Note that you will have to convert this
    to radians.


Once you answer these two question, what you
want to do is to setup the rotation scales to 1.0
and the translation scales to:

      ExpectedRotation / ExpectedTranslation


Note also that, if you need a rotation larger than
(about) 30 degrees, it is unlikely that a registration
will be able to solve it.  You will need to initialize
the Versor with that rotation.


Finally,
you may find useful to look at the following
Insight Journal paper:

http://www.insight-journal.org/browse/publication/304
http://hdl.handle.net/1926/1497
"A Novel Information-Theoretic Point-Set Measure Based on the
Jensen-Havrda-Charvat-Tsallis Divergence"
by Tustison N., Awate S., Gee J.



     Regards,


           Luis



--------------------------------------------------------------------------------------
On Wed, Sep 2, 2009 at 7:21 AM, Neuhaus Jochen <j.neuhaus at dkfz-heidelberg.de
> wrote:

> Hi List,
>
> I need to register two sets of points with unknown point
> correspondences. I use the code of example 2 in the software guide in
> section 8.17, pages 492ff
> (Examples/Patented/IterativeClosestPoint2.cxx), I just exchanged the
> Euler3DTransform with a VersorRigid3DTransform.
>
> When using artificial point sets that have just a uniform translation
> between them, the optimizer finds the correct transformation.
>
> When using real world point sets (from a tracking device) that are also
> rotated and do not fit perfectly, the optimizer fails to find a good
> transformation.
>
> For testing, I established point correspondences manually and used the
> itk::LandmarkTransformInitializer to initialize a
> VersorRigid3DTransform.
> The transform parameters were [0.422982, 0.403292, -0.581769, -27.3493,
> 155.878, 1790.09] (first the 3 versors, then the 3 translation
> parameters). This resulted in a FRE (fiducial registration error) of
> about 1.5mm.
>
> Feeding the same point sets (with correct point correspondences) into
> the optimizer resulted in these transform parameters:
> [2.08919e-006, 0.00995692, 0.116826, -19.5648, 176.277, 1790.94]
>
> As you can see, the versors were hardly changed. The
> EuclideanDistancePointMetric had the following values:
> [56.2236, 51.0406, 28.7475, 109.227, 109.721, 78.9683] (this lead to a
> FRE of about 50mm!)
>
>
> Why does the optimizer stop when the cost function is that bad? As
> proven with the itk::LandmarkTransformInitializer, a much better
> solution exists.
>
> I suspect that the parameter scaling could be the cause. Does anyone
> have a hint of how to set the scales for a VersorRigid3DTransform?
>
>
> My registration code is as following:
>
> /* lots of type definitions */
>  typedef itk::PointSet<mitk::ScalarType, 3> PointSetType;
>
>  typedef itk::EuclideanDistancePointMetric< PointSetType, PointSetType>
> MetricType;
>  typedef itk::VersorRigid3DTransform< double > TransformType;
>  typedef TransformType ParametersType;
>  typedef itk::PointSetToPointSetRegistrationMethod< PointSetType,
> PointSetType > RegistrationType;
>
>  MetricType::Pointer metric = MetricType::New();
>
>  TransformType::Pointer transform = TransformType::New();
>  transform->SetIdentity();
>
>  itk::LevenbergMarquardtOptimizer::Pointer optimizer =
> itk::LevenbergMarquardtOptimizer::New();
>  optimizer->SetUseCostFunctionGradient(false);
>
>  RegistrationType::Pointer registration = RegistrationType::New();
>  // Scale the translation components of the Transform in the Optimizer
>  itk::LevenbergMarquardtOptimizer::ScalesType
> scales(transform->GetNumberOfParameters());
>  const double translationScale = 5000;
>  const double rotationScale = 1.0; // dynamic range of rotations
>  scales[0] = 1.0 / rotationScale;
>  scales[1] = 1.0 / rotationScale;
>  scales[2] = 1.0 / rotationScale;
>  scales[3] = 1.0 / translationScale;
>  scales[4] = 1.0 / translationScale;
>  scales[5] = 1.0 / translationScale;
>  //scales.Fill(0.01);
>  unsigned long numberOfIterations = 80000;
>  double gradientTolerance = 1e-10; // convergence criterion
>  double valueTolerance = 1e-10; // convergence criterion
>  double epsilonFunction = 1e-10; // convergence criterion
>  optimizer->SetScales( scales );
>  optimizer->SetNumberOfIterations( numberOfIterations );
>  optimizer->SetValueTolerance( valueTolerance );
>  optimizer->SetGradientTolerance( gradientTolerance );
>  optimizer->SetEpsilonFunction( epsilonFunction );
>
>  registration->SetInitialTransformParameters(
> transform->GetParameters() );
>    registration->SetMetric( metric );
>  registration->SetOptimizer( optimizer );
>  registration->SetTransform( transform );
>  registration->SetFixedPointSet( targetPointSet );
>  registration->SetMovingPointSet( sourcePointSet );
>
>  registration->Update();
>
>  std::cout << "ICP successful: Solution = " <<
> transform->GetParameters() << std::endl;
>
>
> Sorry for the long email.
> Thank you for your help!
>  Jochen
> _____________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.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
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20090909/99f81173/attachment.htm>