[Insight-users]
Error with MutualInformationHistogramImageToImageMetric
Luis Ibanez
luis.ibanez at kitware.com
Thu Jun 21 18:54:17 EDT 2007
Hi Rosario,
1) Yes, if the error happens before the first optimizer
iteration is reported, then you can safely assume that
the poor initialization of the transform is the source
of the problem.
2) Please look at the example ImageRegistration15.cxx in
the directory: Insight/Examples/Registration
Typically you should compute the translation scaling as
translationScaling[0] = 1 / ( 100 * spacing[0] * size[0] )
translationScaling[1] = 1 / ( 100 * spacing[1] * size[1] )
translationScaling[2] = 1 / ( 100 * spacing[2] * size[2] )
3) You can use Moments if your images are such that if
intensity was equivalent to "mass density", then the
center of mass of both images will be placed in the
same anatomical location. This condition if often not
satisfied in multi-modality registration. In those cases
you should use the Geometry mode, or use the Landmark
based initializer.
4) In principle it will be equivalent to play by multiplying
the step length by a factor K and dividing the parameter
scaling values by a factor of K. Just make sure that when
you consider the product, the values of the versor are
in the range of 0.01. (e.g. about 1 degree rotation)
5) Since you are initializing the versor with a zero angle,
the axis doesn't really matter. This would only be the
source of the problem if your two images are actually
rotated by a large angle. Is that the case ?
6) The selection of the number of histogram bins should be
driven by actually looking at the histograms of both
images. Starting with a 32 bins is a reasonable value.
Note that the number of bins and the number of samples
are somehow related, since you have to make sure that
you have enough samples for providing a decent population
to the number of bins in the joint histogram. In other
words, if you increase the number of bins, you may have
to increase the number of samples accordingly.
Regards,
Luis
==================
Rosario Sance wrote:
>
> Hi Luis,
> Thank you very much for your useful comments. Actually I hadn't
> understood very well the role of the 'maximum step' parameter, and of
> course I have already tried this correction. However, it looked not to
> be enough, since *the error persists*...
> Answering to your questions:
>
> 1) None. Actually I checked that it failed just when beginning the
> registration process. That's why I guess that it was due to
> initialization...
>
> 2) No, but how could I reflect them on the scales? In this particular
> case (but I use to work with different images), the characteristics are:
>
> Size of every volume in the series = 512 x 512 x 24
> ElementSpacing = 0.742188 x 0.742188 x 6.6
>
> 3) I had chosen MomentsOn(), but actually I'd like to know which would
> be more recommendable and which criteria should I consider?
>
> 4) I got it now. But then I have another question. Wouldn't be
> equivalent to change the scales and then allow a larger step? So for
> example (just like an approximation):
>
> /typedef/ / OptimizerType::ScalesType OptimizerScalesType;
> OptimizerScalesType optimizerScales( transform->GetNumberOfParameters() );
> const double translationScale = 1.0 / 1000.0;
> const double rotationScale = 1.0 / 60.0;
>
> optimizerScales[0] = rotationScale ;
> optimizerScales[1] = rotationScale ;
> optimizerScales[2] = rotationScale ;
> optimizerScales[3] = translationScale;
> optimizerScales[4] = translationScale;
> optimizerScales[5] = translationScale;
>
>
>
>
> /I have another 2 comments:
>
> 5) What about the *rotation description?* I didn't change nothing in
> this piece of code, but maybe that's the clue...
>
> /typedef/ / TransformType::VersorType VersorType;
> typedef VersorType::VectorType VectorType;
>
> VersorType rotation;
> VectorType axis;
>
> axis[0] = 0.0;
> axis[1] = 0.0;
> axis[2] = 1.0;
>
> const double angle = 0;
>
> rotation.Set( axis, angle );
> transform->SetRotation( rotation );
>
> /6) Which could be a good size for the *'numberOfHistogramBins' *defined
> for the metric? I am working with unsigned short format, and the images
> have intensities about the range min = 0.000000, max = 660.000000. I
> guess that depending on that parameter the computational cost will
> change, but would you give me a suggestion? 32, 64, 128,... would be any
> of them ok?
>
>
> I'm looking forward to your answer.
> Thanks,
> Rosario
>
>
> At 03:19 13/06/2007, Luis Ibanez wrote:
>
>> Hi Rosario,
>>
>> Thanks for the detailed description of your problem.
>>
>> The combination of registration components that you
>> are using looks fine.
>>
>> You are also guessing correctly that the problem is
>> probably related to the initialization.
>>
>>
>> A couple of questions:
>>
>>
>> 1) How many iterations of the optimizer do you see
>> with the Observer before the Exception is thrown ?
>>
>>
>> 2) I see that in the translation scales you are not
>> taking into account the spacing and number of
>> pixels of your images. It is usually wise to do
>> so because in that way you can be sure that the
>> translation scaling really correspond to the
>> physical extent of your images.
>>
>> What are the number of pixels of your images,
>> and what are their pixel spacings ?
>>
>>
>> 3) What mode of the CenteredTransformInitializer are
>> you using: Geometry ? or Moments ?
>>
>>
>> 4) The value 0.2 as initial maximum step is a bit
>> large. Note that since the parameter scaling
>> for the Versor (unit quaternion) components of
>> the transform is 1.0, that means that your first
>> jump in the parametric space is going to be of
>> 11 degrees in rotation. You may want to start
>> with a smaller amount, maybe something like just
>> 1 degree or rotation, which should be a step of
>> 0.02.
>>
>>
>>
>> Please let us know,
>>
>>
>> Thanks
>>
>>
>> Luis
>>
>>
>>
>>
>> =====================
>> Rosario Sance wrote:
>>
>>> Hi all,
>>> My question is regarding to a particular implementation of a
>>> Registration algorithm I have just tried without success.
>>> Firstly I checked the example " ImageRegistration8.cxx", in which the
>>> use of the " itkVersorRigid3DTransform" is illustrated, and it worked.
>>> Then I did the same with the example "
>>> ImageRegistrationHistogramPlotter.cxx", which deals with the "
>>> MutualInformationHistogramImageToImageMetric", and also worked
>>> properly with my images. (Since they are 3D volumes I did the
>>> required changes and I had no problems at all).
>>> However, when I tried to create a new algorithm with the elements:
>>> TRANSFORM: itkVersorRigid3DTransform (with the
>>> itkCenteredTransformInitializer)
>>> OPTIMIZER: itkVersorRigid3DTransformOptimizer
>>> METRIC: MutualInformationHistogramImageToImageMetric
>>> INTERPOLATOR: itkLinearInterpolateImageFunction
>>> The code was properly built but I achieved a permanent error in
>>> execution:
>>> /ExceptionObject caught !
>>> itk::ExceptionObject (013CF7FC)
>>> Location: "void __thiscall itk::HistogramImageToImageMetric<class
>>> itk::Image<float,3>,class itk::Image<float,3>
>>> >::ComputeHistogram(const class itk::Array<double> &,class
>>> itk::Statistics::Histogram<double,2,class
>>> itk::Statistics::DenseFrequencyContainer> &) const"
>>> File:
>>> c:\itk\insighttoolkit-3.2.0\code\algorithms\itkHistogramImageToImageMetric.txx
>>> Line: 303
>>> Description: itk::ERROR:
>>> MutualInformationHistogramImageToImageMetric(018A50D8): All the
>>> points mapped to outside of the moving image
>>> /I guessed that the problem came from the initialization of the
>>> Optimizer:
>>> /typedef/ / OptimizerType::ScalesType OptimizerScalesType;
>>> OptimizerScalesType optimizerScales(
>>> transform->GetNumberOfParameters() );
>>> const double translationScale = 1.0 / 1000.0;
>>> optimizerScales[0] = 1.0;
>>> optimizerScales[1] = 1.0;
>>> optimizerScales[2] = 1.0;
>>> optimizerScales[3] = translationScale;
>>> optimizerScales[4] = translationScale;
>>> optimizerScales[5] = translationScale;
>>> optimizer->SetScales( optimizerScales );
>>> optimizer->SetMaximumStepLength( 0.2000 );
>>> optimizer->SetMinimumStepLength( 0.0001 );
>>> optimizer->SetNumberOfIterations( 200 );
>>>
>>> /but I was playing with it and I couldn't avoid the fail... The same
>>> with the initialization of the transform:
>>> / / /typedef itk::CenteredTransformInitializer< TransformType,
>>> FixedImageType,
>>> MovingImageType
>>> >
>>> TransformInitializerType;
>>> TransformInitializerType::Pointer initializer =
>>>
>>> TransformInitializerType::New();
>>> (/...)
>>> I have tried it with different images (different sizes and
>>> resolution). All of them are MRI dynamic studies of the abdominal
>>> region (3D volumes on time). I'd really like to know what's happening
>>> with my piece of code. *Is there any specific limitation for mixing
>>> those elements in a new algorithm?
>>> *Thank you very much in advance.
>>> Rosario
>>> ___________________________________________________
>>> Rosario Sance Garzón
>>> Grupo de Tecnologías de Imagen Biomédica
>>> Dpto. Ingeniería Electrónica
>>> E.T.S.I. Telecomunicación - UPM
>>> Ciudad Universitaria
>>> 28040 - Madrid
>>> Teléfono: 91 549 57 00 (Ext. 4220) Fax: 91 336 73 23
>>> Web: http://www.die.upm.es/im
>>> ___________________________________________________
>>> /Rosario Sance Garzón
>>> Biomedical Image Techologies
>>> Electronic Engineering Department
>>> E.T.S.I. Telecomunicación - UPM
>>> Ciudad Universitaria
>>> E-28040 - Madrid (SPAIN)
>>> //Telephone: +34 91 549 57 00 (Ext. 4220) Fax: +34 91 336 73 23
>>> Web: http://www.die.upm.es/im
>>> / ___________________________________________________
>>>
>>> ------------------------------------------------------------------------
>>> _______________________________________________
>>> Insight-users mailing list
>>> Insight-users at itk.org
>>> http://www.itk.org/mailman/listinfo/insight-users
>>
>>
>>
>> --
>> No virus found in this incoming message.
>> Checked by AVG Free Edition. Version: 7.5.472 / Virus Database:
>> 269.8.15/847 - Release Date: 12/06/2007 21:42
>
> ___________________________________________________
>
> Rosario Sance Garzón
> Grupo de Tecnologías de Imagen Biomédica
> Dpto. Ingeniería Electrónica
> E.T.S.I. Telecomunicación - UPM
> Ciudad Universitaria
> 28040 - Madrid
>
> Teléfono: 91 549 57 00 (Ext. 4220)
> Fax: 91 336 73 23
> Web: http://www.die.upm.es/im
> ___________________________________________________
>
> /Rosario Sance Garzón
> Biomedical Image Techologies
> Electronic Engineering Department
> E.T.S.I. Telecomunicación - UPM
> Ciudad Universitaria
> E-28040 - Madrid (SPAIN)
>
> //Telephone: +34 91 549 57 00 (Ext. 4220)
> Fax: +34 91 336 73 23
> Web: http://www.die.upm.es/im
> / ___________________________________________________
>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> 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