[Insight-users] fast marching behavior

Luca Antiga luca.antiga at gmail.com
Tue Jan 12 06:07:00 EST 2010 

Hi Bill,
  indeed this is a "logical" bug, a discrepancy with the way the  
filter implements the original algorithm,
which leads to an output that is not what we expect by reading the  
algorithm.
However, since this class is quite old and its expected use is  
widespread, it's not unreasonable to
think that some of our customers might rely on its uncorrect behavior.  
In many situations the bug
will not affect the solution to a stage where the results are  
incorrect, but they will probably be
numerically different to some extent (for other applications, such as  
Siqi's, things may be worse).
Still, such numerical differences in the solution would lead eventual  
tests or, worse, assumptions based
on expected values in the solution, to fail, and in a very subtle way.
For this reason, I'm all for fixing the bug, but probably a CMake flag  
to allow users to revert to the old
way (that might go away in a couple of releases) is in order, so to  
allow customers a smooth transition.
Or I'm being too conservative here?

Luca



On Jan 12, 2010, at 2:26 AM, Bill Lorensen wrote:

> Luca,
>
> If this is a bug it should be fixed. I don't think we need to maintain
> the bad behavior as long as the API remains the same. Unless I'm
> missing something subtle here.
>
> Perhaps others wish to comment.
>
> Bill
>
> On Fri, Jan 8, 2010 at 11:44 AM, Luca Antiga <luca.antiga at gmail.com>  
> wrote:
>> Hi Dan,
>>  not at all, I was just referring to addressing the latest comments  
>> by Siqi.
>> I have the file on my desktop, waiting for me to handle it :-)
>>
>> Luca
>>
>>
>> On Jan 8, 2010, at 5:18 PM, Dan Mueller wrote:
>>
>>> Hi Luca,
>>>
>>> Is there anything wrong with the "FastMarching.patch" file  
>>> attached to
>>> the bug entry?
>>>   http://public.kitware.com/Bug/view.php?id=8990
>>>
>>> 2010/1/8 Luca Antiga <luca.antiga at gmail.com>:
>>>>
>>>> Hi Siqi,
>>>>  I understand your points. We have to find a solution that  
>>>> guarantees
>>>> backwards compatibility, we cannot switch to setting Alive points
>>>> altogether.
>>>> My proposal is to use the distinction between Trial points and
>>>> InitialTrial
>>>> points in Dan's patch, and to add a flag that will decide whether
>>>> InitialTrial points should be updated or not. This will solve the  
>>>> problem
>>>> maintaining backwards compatibility. With no updating allowed,
>>>> InitialTrial
>>>> points will behave just like you wish, and we won't have to write  
>>>> extra
>>>> code
>>>> for computing the initial layer around Alive points.
>>>> Let me know if you agree with this.
>>>> Dan, I appreciate the lack of time, and I must say we're pretty  
>>>> much on
>>>> the
>>>> same boat. However, this issue carries some weight. As soon as we  
>>>> devise
>>>> an
>>>> acceptable solution, I'm willing to put in some time to commit  
>>>> the patch
>>>> and
>>>> take care of the rest.
>>>> Best regards
>>>> Luca
>>>>
>>>>
>>>>
>>>> On Jan 8, 2010, at 4:06 PM, siqi chen wrote:
>>>>
>>>> The reason I suggest to let the algorithm compute the initial  
>>>> trial is
>>>> the
>>>> following:
>>>>
>>>> If we want to compute a distance map to a perfect circle. The  
>>>> first thing
>>>> we
>>>> do is to discretize the circle, of course, most of the points  
>>>> will be
>>>> non-integer coordinates. To initialize the algorithm, we need to  
>>>> find the
>>>> 4
>>>> neighbor lattice of the discretized circle points and calculate  
>>>> the exact
>>>> distance from these lattice to the circle. In my humble opinion,  
>>>> in the
>>>> current ITK implementation, the algorithm accepts these lattice  
>>>> points as
>>>> "TrialPoints". As I illustrated before, since the FMM method will  
>>>> update
>>>> the
>>>> values of trial points if the new value is smaller, it is very  
>>>> likely
>>>> that
>>>> after the FMM sweeping, the value of the these lattice (the  
>>>> immediate 4
>>>> neighbors of the circle) will change. This will introduce large  
>>>> errors if
>>>> we
>>>> compare the zero iso curve of the result distance to the initial  
>>>> circle.
>>>>
>>>> To fix this problem, there are two options:
>>>> 1. Distinguish between the user input trial and the algorithm  
>>>> generated
>>>> trial as you guys mentioned before. If it's a user input trial,  
>>>> then
>>>> don't
>>>> update the value at all.
>>>> 2. Set these immediate 4 neighbors as KnownPoints, and let the user
>>>> compute
>>>> another layer outside these knownpoints and set them as trial  
>>>> points.
>>>> This
>>>> is quite tedious for the user, since they need to figure out the  
>>>> "upwind"
>>>> themselves. Again, according to the FMM algorithm itself, this  
>>>> initial
>>>> trial
>>>> computation should be done by the algorithm, not the user. That's  
>>>> why I
>>>> think we should let the algorithm to compute the initial trial.
>>>>
>>>> Thanks
>>>> Siqi
>>>>
>>>> On Thu, Jan 7, 2010 at 6:07 PM, Luca Antiga <luca.antiga at gmail.com>
>>>> wrote:
>>>>>
>>>>> Hi Siqi,
>>>>>  this sounds like a perfect contribution for the Insight Journal.
>>>>> I suggest to write a new class for the second-order  
>>>>> implementation.
>>>>> As for the minheap, I agree with you, it would save some memory.
>>>>> As for setting trial points: in the current implementation, if  
>>>>> you don't
>>>>> set trial points you won't have anything to process in the  
>>>>> TrialHeap.
>>>>> I understand your point, but I am kind of in favor of the  
>>>>> current way of
>>>>> working. Setting alive points doesn't necessarily mean that you  
>>>>> want
>>>>> the solution to propagate automatically starting from their  
>>>>> boundary
>>>>> with far points.
>>>>> This could indeed be an option, but it's not
>>>>> necessarily always practically
>>>>> convenient. Is there any particular reason why you wouldn't want  
>>>>> to set
>>>>> your
>>>>> initial points as trial?
>>>>>
>>>>> Luca
>>>>>
>>>>> On Jan 7, 2010, at 11:47 PM, siqi chen wrote:
>>>>>
>>>>> It definitely should be committed. I also have several  
>>>>> suggestions about
>>>>> fastmarching implementation in ITK.
>>>>>
>>>>> 1. I think it's quite misleading to say that "In order for the  
>>>>> filter to
>>>>> evolve, at least some trial points must be specified." I think  
>>>>> it should
>>>>> be
>>>>> "at least some alive points must be specified". The initial  
>>>>> trial points
>>>>> should be calculated by the algorithm itself, not by the user. For
>>>>> example,
>>>>> when I start with FastMarchingImageFilter, I want to try a simple
>>>>> example of
>>>>> calculating a distance map to [50,50], Mr. Kevin Hobbs told me  
>>>>> the right
>>>>> way
>>>>> is to set [50,50] as trial point. However, I think it is not right
>>>>> description. The right way is to set [50,50] as alive points,  
>>>>> and the
>>>>> algorithm computes the initial trial points and start the  
>>>>> iteration.
>>>>>
>>>>> 2. Another thing is about heap implementation. Currently,
>>>>> std::priority_queue is used, and there is some issue as Doxygen
>>>>> described.
>>>>> I think we can use a minheap as the basic data structure. I  
>>>>> wrote a VTK
>>>>> fast
>>>>> maching filter yesterday and I used this min heap data structure.
>>>>>
>>>>> 3. Possible improvements. Current implementation is based on Prof.
>>>>> Sethian' s work, which is actually a first order approximation of
>>>>> distance
>>>>> map. This implementation will introduce large errors in the  
>>>>> diagonal
>>>>> direction. With a little bit change of code, we can implement a  
>>>>> higher
>>>>> order
>>>>> accuracy FMM method. The detail can be found here:
>>>>> http://www2.imm.dtu.dk/pubdb/views/publication_details.php? 
>>>>> id=841 .
>>>>>
>>>>> Thanks
>>>>> Siqi
>>
>> _____________________________________
>> 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
>>

More information about the Insight-users mailing list