[Insight-developers] Bug in ImageFunction causing incorrect results from ResampleImageFilter
Stefan Klein
s.klein at erasmusmc.nl
Mon Jun 23 06:29:02 EDT 2008
Hi Dan,
I'll be happy to work on that!
Let me know once you created the branch.
Kind regards,
Stefan.
Blezek, Daniel J., Ph.D. wrote:
> Stefan,
>
> At the last ITK tcon (just Jim Miller and I), we decided to have a go
> and making points and indices uniform across ITK. Simon Warfield and I
> will be making a branch in a month or so when travel schedules permit.
> The goal is to produce green dashboards with a consistent view of points
> and indices. This will require crawling through much of the
> registration framework, and adapting each ImageFunction to understand
> and appreciate the conventions.
>
> I would welcome any developers interested in this issue.
>
> Regards,
> -dan
>
> -----Original Message-----
> From: Stefan Klein [mailto:s.klein at erasmusmc.nl]
> Sent: Tuesday, June 03, 2008 12:30 PM
> To: Blezek, Daniel J., Ph.D.
> Cc: insight-developers at itk.org
> Subject: Re: [Insight-developers] Bug in ImageFunction causing incorrect
> results from ResampleImageFilter
>
> Hi Dan,
>
> I see your point that the results are counterintuitive if N=1.
>
> But I still think it is a bad idea to just assume some boundary
> condition and increase the valid region. I can imagine that for many
> problems, people would prefer [0,N-1] as a valid range (no need for
> boundary conditions, so simpler/faster code, and no need for any
> arbitrary boundary conditions). By inheriting, more advanced
> interpolators could be implemented, which allow for specifying a
> boundary condition.
>
> Actually, the IsInsideBuffer function definition is always specific for
> an ImageFunction. Only the specific type of ImageFunction (linear
> interpolator, nearest neighbour interpolator, etc) itself knows the
> region on which it can generate a meaningful function value. So, maybe
> it's impossible to define one single implementation that always is
> appropriate.
>
> Regards,
> Stefan.
>
> Blezek, Daniel J., Ph.D. wrote:
>> Hi Stefan,
>>
>> You are indeed correct, the range is closed: [0, N-1]. But I
>> disagree that this is correct behavior. Consider a single slice. It
>> is difficult to have an index inside it because the valid range is: N
>> = 1, [0, N-1] == [0, 0]. Even if my slices are 5mm thick, I am
>> treating them as infinitely thin. In my opinion, this is incorrect.
>> I should reasonably expect to resample a 5mm slice into 5 1mm slices.
>
>> I may not like the results depending on the boundary conditions, but I
>
>> should get "something" out of the resampler. Currently, I would get
>> the first slice, and 4 blank slices.
>>
>> Having our slices go from -0.5 to N-0.5 would match DICOM's
>> definition of a slice. This would be very sensible.
>>
>> OK, forgive the ASCII art... For a two slice volume, interpolated to
>> 10 slices
>>
>> Z 0-----1-----2-----3
>> ITK [+....00000]
>> DICOM [..+....+..]
>> Dan [+....+....]
>>
>> "." are interpolated slices, "+" are slice centers, "0" is where ITK
>> decided the slice was out side the buffer. ITK essentially only gives
>
>> me 5 slices (of data) when I should have 10.
>>
>> 99% of the time we don't care about this issue. When it shows up
>> (like in registration and resampling), it's ugly and we really care.
>>
>> Another issue is backwards compatibility. Do we dare not change
>> this fundamental behavior for fear of breaking existing code? If
>> indeed the community decides to change ITK's definitions, how do we do
>
>> it? Perhaps Bill Lorensen would care to comment?
>>
>> Regards,
>> -dan
>>
>> -----Original Message-----
>> From: Stefan Klein [mailto:s.klein at erasmusmc.nl]
>> Sent: Monday, June 02, 2008 11:39 AM
>> To: Blezek, Daniel J., Ph.D.
>> Cc: insight-developers at itk.org
>> Subject: Re: [Insight-developers] Bug in ImageFunction causing
>> incorrect results from ResampleImageFilter
>>
>> Hi Dan,
>>
>> I'm not sure, but: it seems that with the current implementation the
>> valid range is:
>>
>> [0, N-1]
>>
>> instead of [0,N-1) as you said.
>>
>> (if index==m_EndContinuousIndex the method returns true).
>>
>> In my opinion, the current implementation is right. For example, I see
>
>> no meaningful default implementation for a linear interpolator at
>> indices > N-1. (unless some boundary condition has been set).
>>
>> Moreover, related to the discussion on the voxel centre (IndexToPoint
>> conversion): wouldn't it be more logical to define the valid region as
>
>> going from -0.5 to N-0.5?
>> (but that would cause troubles again for some interpolators, so i
>> would vote to keep it like it is right now)
>>
>> Kind regards!
>> Stefan
>>
>>
>> Blezek, Daniel J., Ph.D. wrote:
>>> I have some fairly thick CT slices (5mm) that I would like to
>>> resample
>>> to be isotropic. A straightforward use of ResampleImageFilter leaves
>>> 4 blank slices in the output volume. I've tracked this down to a
>>> bug(?) in ImageFunction. This code is incorrect IMHO:
>>>
>>> virtual bool IsInsideBuffer( const ContinuousIndexType & index )
>> const
>>> {
>>> for ( unsigned int j = 0; j < ImageDimension; j++ )
>>> {
>>> if ( index[j] < m_StartContinuousIndex[j] ) { return false;
> };
>>> if ( index[j] > m_EndContinuousIndex[j] ) { return false; };
>>> }
>>> return true;
>>> }
>>>
>>> (Dealing with Z in 3D) This code assumes the last slice has zero
>>> thickness. For N slices, the range where IsInsideBuffer returns true
>
>>> is [0, N-1). This is incorrect, the correct result should be [0, N).
>>> [0,
>>> N) gives the last slice it's proper thickness. The windowed Sinc
>>> interpolator properly handles the last slice, i.e. [N-1, N), but I
>>> have not checked others. A proper implementation would look like
>>> this
>>> (note the vcl_floor calls):
>>>
>>> if ( vcl_floor(index[j]) < m_StartContinuousIndex[j] ) {
>>> return false; };
>>> if ( vcl_floor(index[j]) > m_EndContinuousIndex[j] ) {
>>> return false; };
>>>
>>>
>>> I've run all the tests on a fresh check out with the above change.
>>> As
>>> may be expected, almost all of the registration tests fail, while
>>> most
>>> of the resample tests pass. Surprisingly,
>>> ResampleVolumesToBeIsotropicTest passed...
>>>
>>> Any comments on this bug? Should the tests be made to accommodate
>>> the
>>> change, or should I override the method in a custom interpolator?
>>>
>>> Regards,
>>> -dan
>>>
>>>
>>>
>>> Daniel Blezek, PhD
>>> Medical Imaging Informatics Innovation Center Enterprise Imaging
>>> Systems Unit
>>>
>>> P 127 or (77) 8 8886
>>> T 507 538 8886
>>> E blezek.daniel at mayo.edu
>>>
>>> Mayo Clinic
>>> 200 First St. S.W.
>>> Harwick SL-44
>>> Rochester, MN 55905
>>> mayoclinic.org
>>>
>>> _______________________________________________
>>> Insight-developers mailing list
>>> Insight-developers at itk.org
>>> http://www.itk.org/mailman/listinfo/insight-developers
More information about the Insight-developers
mailing list