[ITK] [ITK-users] Inverse of Versor Rigid Transform and RIRE
Matt McCormick
matt.mccormick at kitware.com
Mon Apr 13 12:38:56 EDT 2015
Hi Gabriel,
Good to hear about your progress!
To get a better idea of what is happening in the optimization, it may
insightful to visualize the cost function and where the optimizer is
moving over that space. There may be that there is a local minima that
is catching the optimizer. See, for example, the later parts of this
video, which displays the cost function as a surface and the
progression of the optimization [1]. See the Metrics section of the
ITK Software Guide for another example and discussion on how to do
this. We can make use of the ITK Command class, the IterationEvent,
and the ExhaustiveOptimizer.
The code itself is the best place to find details on multithreading in
the metrics. As a first step, you can call SetNumberOfThreads on an
Object or call itk::MultiThreader::SetGlobalDefaultNumberOfThreads.
There is a lot of interest in improving the mutual information metric
performance -- please keep the list informed of your experiences.
Thanks,
Matt
[1] https://drive.google.com/file/d/0B986LSX8iqF-bl9qQWZobmRWeEk/view?usp=sharing
On Thu, Apr 9, 2015 at 9:17 PM, Gabriel A. Giménez
<gabrielgimenez85 at gmail.com> wrote:
> Hi Matt!
>
>> If you are working with a current version of ITK, when you get the
>> inverse transform, it should transfer the Center for you. The Center
>> point locations are the FixedParameters for the parent class,
>> MatrixOffsetTransformBase. The relationship with parent classes can be
>> found by examining the Doxygen page for the class [1].
>
>
> I took a few days to read the RIRE documentation and I came to the
> conclusion that the direction of registration is not really
> relevant...taking up CT as the fexed image and MR as the moving image can be
> used the transformation that provides ITK...without having to calculate the
> inverse.
>
>
>> Do you get the same result by applying the TransformPoint() method of
>> the inverse transform? This is the API call that should be applied to
>> transform a point.
>
>
> Yes Matt, I get the same results using the API.
>
> But...I made some changes in my optimizers and I got very good results, like
> this ( with RegularStepGradientDescentOptimizerv4 and
> MattesMutualInformationImageToImageMetricv4):
>
> iterations = 200
> Metric value = -0.709388
> versor X = 0.0155583
> versor Y = 0.00642035
> versor Z = -0.0487144
> Translation X = 7.82977
> Translation Y = -60.1034
> Translation Z = -23.6258
> +-----------------------------------------------------------------------------------------------+
> | X GT| Y GT| Z GT| X R|
> Y R| Z R|
> -----------------------------------------------------------------------------------------------+
> | -7.573100| -41.253400| -27.309300| -7.661395|
> -40.915138| -26.044441|
> | 324.872200| -72.815900| -32.906300| 324.712907|
> -73.345104| -30.833635|
> | 24.160700| 291.039300| -16.272700| 24.902019|
> 291.325008| -15.874605|
> | 356.606000| 259.476800| -21.869700| 357.276322|
> 258.895042| -20.663798|
> | -6.055400| -45.115700| 84.613700| -6.394922|
> -44.465633| 85.892103|
> | 326.389900| -76.678200| 79.016800| 325.979381|
> -76.895599| 81.102910|
> | 25.678400| 287.176900| 95.650300| 26.168493|
> 287.774513| 96.061940|
> | 358.123700| 255.614500| 90.053400| 358.542796|
> 255.344547| 91.272747|
> +-----------------------------------------------------------------------------------------------+
> [X, Y, Z]GT are the "ground truth" values and [X, Y, Z]R are my results
>
> Now, something I find strange is that when increasing the number of
> iterations...metric value limprovement is too little but the result is
> little worse..., example:
>
> Iterations = 334
> Metric value = -0.710918
> versor X = 0.0216566
> versor Y = 0.00700629
> versor Z = -0.0508766
> Translation X = 7.80722
> Translation Y = -60.5124
> Translation Z = -24.1047
>
> +-----------------------------------------------------------------------------------------------+
> | X GT| Y GT| Z GT| X R|
> Y R| Z R|
> +-----------------------------------------------------------------------------------------------+
> | -7.573100| -41.253400| -27.309300| -8.342271|
> -39.764911| -28.121895|
> | 324.872200| -72.815900| -32.906300| 323.882938|
> -73.594962| -33.530625|
> | 24.160700| 291.039300| -16.272700| 25.690487|
> 292.179801| -13.916412|
> | 356.606000| 259.476800| -21.869700| 357.915696|
> 258.349750| -19.325143|
> | -6.055400| -45.115700| 84.613700| -7.022108|
> -44.688304| 83.762051|
> | 326.389900| -76.678200| 79.016800| 325.203101|
> -78.518355| 78.353321|
> | 25.678400| 287.176900| 95.650300| 27.010650|
> 287.256408| 97.967534|
> | 358.123700| 255.614500| 90.053400| 359.235859|
> 253.426357| 92.558803|
> +-----------------------------------------------------------------------------------------------+
>
> This pattern is repeated with other optimizers ( like OnePlusOne and a GA
> approach that I am implementing ), that you think about it?
>
> Other questions Matt...:
>
> How works the multithreaded in metrics ? is customizable? improve
> performance? specifically in the case of Mattes Mutual Information...
>
> I tried using the helper CenteredVersorTransformInitializer... but the
> transformation that generates makes, incredibly and also very strange,the
> optimizers does not advance...using CenteredTransformInitializer this does
> not happen...
>
> Really Thanks in advance Matt!
> Regards,
>
>
>
> 2015-04-07 15:03 GMT-04:00 Matt McCormick <matt.mccormick at kitware.com>:
>
>> Hi Gabriel!
>>
>> > I am use RIRE project, specifically CT (movig) and MR_PD (fixed)
>> > images.
>> > Basically, I hava a set of point (in millimeters) of the CT image to
>> > which
>> > apply the trasform result of the registration and updaload this results
>> > in
>> > the web for the evaluation. Example of set of points and his "ground
>> > truth"
>> > :
>> >
>> > Point x y z new_x new_y new_z
>> >
>> > 1 0.0000 0.0000 0.0000 -7.5731 -41.2534 -27.3093
>> > 2 333.9870 0.0000 0.0000 324.8722 -72.8159 -32.9063
>> > 3 0.0000 333.9870 0.0000 24.1607 291.0393 -16.2727
>> > 4 333.9870 333.9870 0.0000 356.6060 259.4768 -21.8697
>> > 5 0.0000 0.0000 112.0000 -6.0554 -45.1157 84.6137
>> > 6 333.9870 0.0000 112.0000 326.3899 -76.6782 79.0168
>> > 7 0.0000 333.9870 112.0000 25.6784 287.1769 95.6503
>> > 8 333.9870 333.9870 112.0000 358.1237 255.6145 90.0534
>>
>> Trying to reproduce previous results is a good path forward.
>>
>> > So, the first I need is the transformation to apply, for that I do the
>> > following :
>> >
>> > //get the inverse transform
>> > TransformType::Pointer inverseTransform = TransformType::New();
>> > inverseTransform->SetCenter( finalTransform->GetCenter() );
>> > bool response = finalTransform->GetInverse(inverseTransform);
>> >
>> >
>> > It makes sense to use the same center in the inverse transform?. A
>> > "quaternion" define an "axis" (right part) of rotation and an angle to
>> > use
>> > for rotate the image about this axis...why use a center of rotation...?
>>
>> If you are working with a current version of ITK, when you get the
>> inverse transform, it should transfer the Center for you. The Center
>> point locations are the FixedParameters for the parent class,
>> MatrixOffsetTransformBase. The relationship with parent classes can be
>> found by examining the Doxygen page for the class [1].
>>
>>
>> > Second, apply this transform...as follows:
>> >
>> > NewPoint = RotationMatrix * OriginalPoint + Offset
>> >
>> > The rotation matrix and the offset are obtained from the inverse
>> > transforme
>> > objetc. Found something wrong? something that is not taking into account
>> > ?
>> > The results do not appear to be correct...the calculated error is too
>> > big
>> > and does not correspond with the visual result.
>>
>> Do you get the same result by applying the TransformPoint() method of
>> the inverse transform? This is the API call that should be applied to
>> transform a point.
>>
>> Thanks,
>> Matt
>>
>>
>> [1]
>> http://www.itk.org/Doxygen/html/classitk_1_1VersorRigid3DTransform.html
>
>
>
>
> --
> Gabriel Alberto Giménez.
_____________________________________
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.php
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://public.kitware.com/mailman/listinfo/insight-users
More information about the Community
mailing list