[Insight-users] Subtracting Images - A Possible Pixel-Type Issue

[Luis Ibanez]

Hi Michael,

Thanks for the detailed description of the problem.

By looking at the code of the Functor in the
itkSquaredDifferenceImageFilter:

   inline TOutput operator()( const TInput1 & A,
                              const TInput2 & B)
   {
     const double dA = static_cast<double>( A );
     const double dB = static_cast<double>( B );
     const double diff = dA - dB;
     return static_cast<TOutput>( diff * diff );
   }


The only way this can generate negative values
is by overflow in the range of the TOutput=short
type.

You probably have in your images difference that
are large enough to overflow a short when the
values are put to the square. Eg. A short overflows
at "32767", which means that a difference of "181"
(when elevated to the square) is enough to saturate
the output and to roll over to the negative range.

For example, if the difference between the pixel
values of two images is = 183, then its square
will be: 33493, and when stored in a signed short
it will become = 33493-32767 = -721.

------

What you may want to do is to replace this filter
with the itkAbsoluteValueDifferenceImageFilter,
whose Functor does:

   inline TOutput operator()( const TInput1 & A,
                              const TInput2 & B)
   {
     const double dA = static_cast<double>( A );
     const double dB = static_cast<double>( B );
     const double diff = dA - dB;
     const double absdiff = ( diff > 0.0 ) ? diff : -diff;
     return static_cast<TOutput>( absdiff );
   }


and therefore will remain in the range of shorts.

Note that if you are out of luck (not to bring
Murphy's Laws into this discussion... but...)
you can still saturate a short in this context.

If you have two pixels:

    P1 = +20000
    P2 = -20000

which are valid short values,
their subtraction will be computed incorrectly
in short:

    DP = 40000

will be stored as:  40000-32767 = -7232


Which means that if you really care about the
safety of the filter, you should customize your
own with some clamping rule such as:

   inline TOutput operator()( const TInput1 & A,
                              const TInput2 & B)
   {
     const double dA = static_cast<double>( A );
     const double dB = static_cast<double>( B );
     const double diff = dA - dB;
     const double absdiff = ( diff > 0.0 ) ? diff : -diff;
     if( absdiff > itk::NumericTraits<TOutput>::max() )
       {
       absdiff = itk::NumericTraits<TOutput>::max();
       }
     return static_cast<TOutput>( absdiff );
   }


You can easily write this filter by editing the
.h file of either one of the two filters described
above.

Of course, if you write this new filter, it will
be great is you submit it as a contribution to the
Insight Journal:

http://www.insight-journal.org/InsightJournalManager/index.php




Please let us know if you run into any problems,


    Thanks


       Luis


------------------
M. Wirtzfeld wrote:
> 
> Hello all,
> 
> I am relative new to ITK, so bare with me...
> 
> 
> I am working on a registration project and I would like to take the 
> difference between an image slice and a resulting image taken from a 
> volume after registration done.  In essence, a pixel-by-pixel subtraction.
> 
> 
> The itk::SquaredDifferenceImageFilter is currently being used to 
> determine the difference image.  However, when pixel-values are reviewed 
> in the resulting image, large negative values are found.  It is believed 
> these negative values are due to the use of the SHORT data-type to 
> represent the pixel-type.  The squared-operation of the filter obviously 
> handles negative subtraction results;  likely to occur regardless of the 
> order in which the two images are given to this filter.
> 
> Is the use of the SHORT data-type likely responsible for these negative 
> pixel-values?  Changing the type is not an option.
> 
> 
> Reading through the ITK Software Guide (Version 2.4), I came across 
> the itk::SubtractImageFilter topic in Section 8.2, "Hello World" 
> Registration, page 323.  In this section it says,
> 
> 
> "... the use of subtraction as a method for comparing the images is 
> appropriate here because we chose to represent the images using a pixel 
> type float.  A different filter would have been used if the pixel type 
> of the images were any of the unsigned integer type."
> 
> 
> 
> The difference is defined, with reference to respective class 
> documentation, as "Output = Input1 - Input2" for the 
> itk::SubtractImageFilter.  With this in mind, and the fact that the 
> pixel-values across the entire image for "Input2" will likely not be 
> less than the respective pixels in "Input1", can I directly apply the 
> itk::SubtractImageFilter in this case?
> 
> 
> 
> Thank you,
> 
> Michael.
> Robarts Research.
> 
> 
> ------------------------------------------------------------------------
> 
> _______________________________________________
> Insight-users mailing list
> Insight-users at itk.org
> http://www.itk.org/mailman/listinfo/insight-users