[Insight-users] Interpolation and centered pixel coordinates

[Luis Ibanez]

Hi Dženan,

Thanks for sharing your findings.

This is quite interesting.

Did you log a bug entry for this problem ?

I'm trying to find the exact place with this
conversion of:

double cosFi=vec1*vec2/(vec1.length()*vec2.length());

Is this something that you found in ITK ?
or is this part or your application's code ?


   Please let us know,


           Thanks


                Luis


------------------------------------------------------
2010/6/16 Dženan Zukić <dzenanz at gmail.com>

> Hi Luis,
>
> I found out what the problem is.
>
> The way I call evaluate is OK:
>
> image->TransformPointFromPhysicalToContinouosIndex(p, ind, ok);
> if (ok)
>   pixel=interpolator->EvaluateAtContinuousIndex(ind);
>
> but my original point p had NaN values in it, which were caused by code
> which looked something like this:
>
> double cosFi=vec1*vec2/(vec1.length()*vec2.length());
> //cosine of angle between vectors
> p+=acos(cosFi)*...;
>
> Due to numerical round-off errors it would rarely occur that cosFi>1
> (something like 1.00000002534), and acos would produce NaN, which would then
> propagate all the way to poor interpolator which would then complain.
> Shortly after starting to use interpolator I also started using different
> sort of input files, which made my program do much more iterations of the
> above numerically unstable vector code and made manifestations of that
> problem more often.
>
> Since the problem was manifesting itself in interpolator, and just before
> calling interpolator I was getting its input data from
> TransformPhysicalPointToContinuousIndex, and code for this transformation is
> very convoluted, I suspected it might be wrong and I had to ask.
>
> I am also quite busy these last few weeks (I barely managed to find time to
> fix that bug which I adopted :D), so I soon forgot I asked this question,
> and of course didn't write that I found out the problem (which was quite
> far, in calculation of normals and curvatures of the mesh).
>
> Thanks for the answer,
> Dženan
>
> 2010/6/16 Luis Ibanez <luis.ibanez at kitware.com>
>
>
>> Dženan,
>>
>> We did quite a heavy testing when we fixed
>> the problem of centered pixel coordinates.
>>
>> The Dashboard was green when we released
>> ITK 3.18, which give us some confidence on
>> the correctness of the code.
>>
>> --
>>
>> Please note that it is your responsibility to call
>>
>>      interpolator->IsInsideBuffer( point )
>>
>> before you call
>>
>>      interpolator->Evaluate( point );
>>
>>
>> That is, your code should take the form:
>>
>>
>>    if(    interpolator->IsInsideBuffer( p1 ) )
>>       {
>>       p2 = interpolator->Evaluate( p1 );
>>       }
>>
>>
>> --
>>
>> If you still find run-time errors, please
>> post to the list a minimal code example
>> that displays the problem.
>>
>>
>>     Thanks
>>
>>
>>         Luis
>>
>>
>> ----------------------------------------------------------
>> 2010/6/11 Dženan Zukić <dzenanz at gmail.com>
>>
>>> Hi everyone,
>>>
>>> I am running into some crashes in
>>> \InsightToolkit-3.18.0\Code\Common\itkLinearInterpolateImageFunction.txx
>>> (line 150), and I wonder could it be due to recent move to consistent use of
>>> centered pixels (physical coordinates going from -spacing/2 to
>>> (n-1)*spacing+spacing/2, instead of 0 to n*spacing). Did anyone have
>>> encounters with similar problems? More specifically, does
>>> TransformPhysicalPointToContinuousIndex in
>>> \InsightToolkit-3.18.0\Code\Common\itkImageBase.h honor the centered pixel
>>> coordinates?
>>>
>>> Regards,
>>> Dženan
>>>
>>> _____________________________________
>>> 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.html
>>>
>>> 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://www.itk.org/mailman/listinfo/insight-users
>>>
>>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20100617/31c7e354/attachment-0001.htm>