[Insight-users] Not so Silly questions about registration

Luis Ibanez luis.ibanez@kitware.com
Tue, 25 Feb 2003 10:53:18 -0500


Hi Daniel,

I'm affraid to have to disagree.

It doesn't make perfect sense to replace floats with chars
and expect to get the same results, not even with simple
operations.

Let's consider:

     T a;
     T b;
     T c = a * b;
     T d = a / b;
     T e = a - b;
     T f = a + b;

Can we expect this to work the same for T=char and T=float ?

T=float and T=double will already give different results...

and if you use VC++, T=float in debug mode and T=float in
release mode will also give different results.

Even if we use the NumericTraits<> and cast to the computational
types, there are always lateral effects the we have to care about.

Do we want to discover these secondary effects at run time in
an image guided surgery system ?

or do we want the developers to deal with the problem at compile
time when they can afford the time to analyze the problem and fix
it while having enough coffee available.

If we put ourselves in the shoes of a naive user of Insight,
the first thing we should do is to take the naivety out of him.
Naivety is not a desirable quality for a developer of medical
image applications.

Insight is not only a software toolkit, it is also an educational
resource. I think we should be commited to insist in the special
requirements of medical imaging as opposed to plain computer
graphics concepts. Insight is not Photoshop. The images that we
are intended to process contain information about somebody lying
on a bed in a hospital. Taking care of the implications of data
types doesn't seem to be such a burden when we consider it in
this perspective.

--

Unfortunately we cannot restrict the registration classes to work
only on float images. The reason is that some registration problems
can actually be solved in char pixel types.  It depends on the metric
we use, and users have the freedom to define their own metrics.

As in real life, Freedom comes with Responsibility attached.

Insight is a toolkit, there is a limit to the things we can do
for preventing somebody from using a hammer to take out a screw.
We could install an anti-screw security system in the hammer but
this will probably make the hammer useless for hitting nails.
This is, in my opinion, more of a documentation and education
issue than a programming one.

As developers we cannot be charged with getting the same results
across a range of datatypes. That is simply not possible. Can
we add two constant images with values 250 in unsigned char ?
There are, after all, good reasons for different types to exist.

The fact that we can instantiate it, doesn't mean it will work.

We can instantiate an image of PixelType = std::string. No problem
at all, but we cannot expect a median filter to work on this image,
nor to be able to save it in a PNG file.

----

The GradientMagnitudeRecursiveGaussianImageFilter requires so much
memory because it provides a trade-off between memory and speed.
If memory becomes a problem, then you have the option of paying in
speed and use a convolution based filter like the 
GradientMagnitudeImageFilter. It will take 10 times less memory
but it will probably run 10 times slower.

----


It is always healthy to revisit the frameworks and verify if they
can be improved. That's after all the main advantage of an Open
Source project. I will be glad to analyze the options that could
lead to solve registration problems on a whole range of data
types and discuss their possible implementations.

Proabably the best way to start is to create the trivial test you
suggested in which we could evaluate which types behave correctly
for different registration setups.


Would you like to give it a try to implement this test ?
I'll be happy to help.


Best regards



    Luis



---------------------------------------------

Blezek, Daniel J (Research) wrote:
> Luis, Dave,
> 
>   I'm not sure this is the best way of addressing Dave's concerns.  Let's put ourselves in the shoes
> of a naive user of Insight.  Dave did exactly what I would have done.  He said, "Hmm, why waste 4x of
> my memory by using floats, I'll just run on unsigned chars, and get the same result.  This makes
> perfect sense (to me at least).
> 
>   IMHO, it's OK to be PixelType-Aware, but Dave would have spent an amazing amount of time trying to
> figure out why he got different results changing from float to unsigned char for registration.  While
> I agree that floats _can_ be faster, you can also be killed if you start swapping.  If we really want
> to catch problems at compile time, shouldn't the registration classes prevent you from instantiating
> a non-floating point class?  Again, as a naive user, I would think one of the best benefits of ITK is
> being able to create an unsigned char version of MI registration if I want, but apparently, that is
> not the case.  Shouldn't the ITK developer be charged with getting the same results across a range of
> datatypes?  How hard would it be to make the GradientRecursiveGaussianImageFilter have an output of
> float, and do all internal calculations with floats/doubles?
> 
>   I was a bit put off that the GradientRecursiveGaussianImageFilter requires 10x the image size to do
> it's work.  That seems really huge.  For example, if I have a 512x512x512 CT (268M, fixed) to
> register to a 256x256x128 MR image (16M, moving).  I'd end up having to convert to float ((268+16) *
> 2 = 568M), and then the GradientRecursizeGaussianImageFilter would take 10x 32M = 320M, for a grand
> total of 888Meg of ram?  This would send most systems into swap death very quickly.
> 
>   Should the registration hierarchy be re-visited and be made type independant?  It would be trivial
> to make a test that would register two images across a whole range of data types and verify the same
> results are obtained.
> 
> -dan
> 
> --
> Daniel Blezek, Ph.D.
> blezek@crd.ge.com
> Visualization and Computer Vision Lab, Imaging Technologies
> GE Global Research Center
> 
> 
> 
>>-----Original Message-----
>>From: Luis Ibanez [mailto:luis.ibanez@kitware.com]
>>Sent: Monday, February 24, 2003 4:01 PM
>>To: holmesd3@yahoo.com
>>Cc: lng@insightful.com; insight
>>Subject: Re: [Insight-users] Silly questions about registration
>>
>>
>>
>>Hi David
>>
>>
>>Your experience just hightlight the importance of
>>being PixelType-Aware in Medical Image applications.
>>
>>ITK methods are pixel-type templated but this doesn't
>>mean they are type-independent.
>>
>>On the contrary, ITK imposes a higher responsibility
>>on the developer to verify the consequences of his/her
>>data-type decisions.
>>
>>The great advantage is that by doing this type checking
>>at compile type, you will be confronted with any potential
>>problems early in the development cycle, instead of being
>>surprised when you are presenting a final demo of your
>>system.
>>
>>---
>>
>>Using a "unsigned char" image as input will not
>>result in the same output during registration.
>>
>>The reason is that the MeanSquaresImageMetric
>>computes its derivatives by multiplying the Jacobian
>>matrix of the transformation with the gradient of
>>the fixed image.
>>
>>This image gradient is computed with the
>>GradientRecursiveGaussianImageFilter which uses the
>>type of the input image for storing intermediate
>>results. Using an "unsigned char" image as input
>>results in a lower quality computation of the
>>gradient since the intermediate values are
>>systematically truncated.
>>
>>Poor image gradients lead to poor metric derivatives.
>>
>>Poor metric derivatives lead to poor updates of the
>>position in the parametric space by the optimizer.
>>
>>Poor updates in the parametric space result in incorrect
>>values for the final transform. In your case, it seems
>>to simply send the transform in a divergent path.
>>
>>
>>
>>The memory savings that you may achieve with the
>>PixelType replacement have to be evaluated in the
>>context of your application. For example, the
>>GradientRecursiveGaussianImageFilter will take
>>about 10 times the memory of the input image
>>regardless of the input image type.
>>
>>Your application will not run faster for using
>>unsigned char as input. In fact it may even run
>>slower. Float types will be managed by the math
>>instruction of the processor which are faster
>>than math operations on integers. Doing math with
>>unsigned chars imposes an extra burden on the
>>processor and wastes capacity in data transfer
>>in-out the processor.
>>
>>If speed is a concern, you should rather consider
>>taking advantage of the multiresolution framework
>>and the transform initializers. These strategies
>>will provide a more effective way of reducing
>>computational times.
>>
>>
>>Please let us know if you have further questions.
>>
>>
>>  Thanks
>>
>>
>>   Luis
>>
>>
>>-----------------------------------
>>
>>David Holmes wrote:
>>
>>
>>>I am unfortunately still having some trouble with the
>>>registration methods.  In particular, I am still
>>>working with the ImageRegistration1.cxx example.  When
>>>I run the code as is, it gives the correct output and
>>>I am happy.  However, when looking and the code and
>>>data, I realized that the data is 8-bit data but in
>>>the code we are working with floats.  In the interest
>>>of memory management and speed, I decided to change
>>>the typedef from:
>>>
>>>typedef   float   PixelType
>>>
>>>to:
>>>
>>>typedef  unsigned char   PixelType
>>>
>>>Because I left everthing else the same, I assumed that
>>>the algorithm would work simply with the unsigned char
>>>data.  When I ran it, I got:
>>>
>>>Result = 
>>> Translation X = -191.882
>>> Translation Y = -190.493
>>> Iterations    = 70
>>> Metric value  = 13.5914
>>>
>>>I haven't dug into the ITK code nor have I added an
>>>Observer to the optimizer yet.  I will do that as a
>>>next step, but I figured that someone might have a
>>>better idea about why this occurse before I dig more
>>>myself.
>>>
>>>Thanks
>>>
>>>david
>>>
>>
>>
>>_______________________________________________
>>Insight-users mailing list
>>Insight-users@public.kitware.com
>>http://public.kitware.com/mailman/listinfo/insight-users
>>
> 
>