[Insight-users] fast marching behavior
Luca Antiga
luca.antiga at gmail.com
Thu Jan 7 18:07:24 EST 2010
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
>
>
> On Thu, Jan 7, 2010 at 4:56 PM, Luca Antiga <luca.antiga at gmail.com>
> wrote:
> Dan, Siqi,
> the patch should be committed. I remember the original discussion,
> we never followed up
> committing the patch. It's time to do it.
> However, since the fast marching method is probably used by a lot of
> people, how about
> adding a cmake flag to revert to the "wrong" implementation to
> ensure backward compatibility?
>
>
> Luca
>
>
> On Jan 7, 2010, at 10:42 PM, siqi chen wrote:
>
>>
>> I think the patch correctly implements the FMM method described by
>> Sethian.
>>
>> Thanks
>> Siqi
>>
>>
>> On Thu, Jan 7, 2010 at 1:54 PM, Dan Mueller <dan.muel at gmail.com>
>> wrote:
>> Hi Siqi,
>>
>> I managed to find some time to take a look at the issue you reported.
>>
>> I think the following patch will help:
>> http://public.kitware.com/Bug/view.php?id=8990
>>
>> Here is the results I get with the patch applied:
>>
>> The distance from [50,50] to [49.7,49.8] is: 0.36
>> The distance from [49,49] to [49.7,49.8] is: 1.06
>> The distance from [49,50] to [49.7,49.8] is: 0.73
>> The distance from [50,49] to [49.7,49.8] is: 0.85
>> The distance from [48,49] to [49.7,49.8] is: 2.06
>> The distance from [48,50] to [49.7,49.8] is: 1.73
>> The distance from [49,48] to [49.7,49.8] is: 2.06
>> The distance from [49,51] to [49.7,49.8] is: 1.73
>> The distance from [50,48] to [49.7,49.8] is: 1.85
>> The distance from [50,51] to [49.7,49.8] is: 1.36
>> The distance from [51,49] to [49.7,49.8] is: 1.85
>> The distance from [51,50] to [49.7,49.8] is: 1.36
>>
>> Is that what you expected?
>>
>> Hope this helps.
>>
>> Cheers, Dan
>>
>> 2010/1/5 siqi chen <siqichensc at gmail.com>:
>> >
>> > I read a few more papers and thought about the question number 2
>> again. I
>> > think there is misinterpretation on my part. The 4 neighbors of
>> [49.7,49.8]
>> > along their exact distance value should be set as Alive Points
>> instead of
>> > Trial Points. Then the neighbors of the 4 neighbors (another
>> layer around
>> > the 4 neighbors) are set as trial points, their initial tentative
>> values are
>> > calculated using upwind difference method. When I check the
>> result distance
>> > map, a few things to notice:
>> > 1. The distance value of the 4 neighbors of [49.7, 49.8] are
>> correct. This
>> > is obvious since I set them as Known instead of Trial.
>> > 2. The distance value of the neighbors of the 4 neighbors (the
>> input trial
>> > points) still changed. I think this is due to the inherent FMM
>> accuracy and
>> > update method.
>> >
>> > You can try this with the following code.
>> > http://www.rpi.edu/~chens/download/main3.cpp
>> >
>> > Thanks
>> > Siqi
>> >
>> > On Mon, Jan 4, 2010 at 6:03 PM, siqi chen <siqichensc at gmail.com>
>> wrote:
>> >>
>> >> To better illustrate my questions regarding fast marching, I put
>> 2 example
>> >> code in the attachment.
>> >>
>> >> In main1.cpp, I simply compute a distance map to point [50,50].
>> As you can
>> >> see from the output, the distances from the 4 neighbors of
>> [50,50] to
>> >> [50,50] are correct, obviously the result is 1. However, the
>> distance from
>> >> [51,51] to [50,50] is 1.707 instead of 1.414, which is obviously
>> wrong. I
>> >> think this is due to the fast marching accuracy itself. If we
>> switch to
>> >> higher order FMM, the result should be improved.
>> >>
>> >> In main2.cpp, I perturb the target point a little bit. Instead,
>> I want to
>> >> compute the distance map to point [49.7,49.8]. From my point of
>> >> understanding, I need to initialize the 4 neighbors of [49.7,
>> 49.8] and put
>> >> them into the TrialPoints. As you can see, I compute the exact
>> distance from
>> >> these 4 neighbors to [49.7,49.8] and put them into TrialPoints.
>> However,
>> >> when I go back and check the result distance map, some thing is
>> different.
>> >> The distances from these 4 neighbors to [49.7,49.8] are changed.
>> As you can
>> >> see, the distance from [50,50] to [49.7,49.8] remains correct.
>> This is
>> >> because this value is the smallest in the TrialPoints, therefore
>> it is
>> >> pushed in to the AlivePoints heap first and the value is frozen
>> since then.
>> >> I think there is something wrong here about whether to update
>> trial points
>> >> value or not. If this trial point is user specified, then the
>> value should
>> >> not be updated. I noticed a related discussion a couple of
>> months ago in the
>> >> mailing list,
>> >> http://www.itk.org/pipermail/insight-users/2009-May/030282.html
>> >>
>> >> http://www.rpi.edu/~chens/download/main1.cpp
>> >> http://www.rpi.edu/~chens/download/main2.cpp
>> >>
>> >> Any input is appreciated.
>> >> Siqi
>> >>
>> >>
>> >> On Mon, Jan 4, 2010 at 3:32 PM, Dan Mueller <dan.muel at gmail.com>
>> wrote:
>> >>>
>> >>> Hi Siqi,
>> >>>
>> >>> Indeed I am familiar with Fast Marching. I saw your question to
>> the
>> >>> mailing list, but did not respond because I have not
>> experienced what
>> >>> you describe: when I set the trial point value, that is the
>> value in
>> >>> the arrival function.
>> >>>
>> >>> Perhaps you could post to the mailing list a minimal example
>> >>> (code+cmake+data) demonstrating your issue. That would make it
>> really
>> >>> easy for me to help you!
>> >>>
>> >>> Cheers, Dan
>> >>>
>> >>> 2010/1/4 siqi chen <siqichensc at gmail.com>:
>> >>> > Hi, Dan,
>> >>> >
>> >>> > Sorry to bother you. From the ITK mailing list, I noticed you
>> reported
>> >>> > a bug
>> >>> > about FastMarchingImageFilter couple of months ago. So I
>> guess you are
>> >>> > a
>> >>> > fast marching expert : )
>> >>> >
>> >>> > I am trying to use FastMarchingImageFilter to calculate a
>> distance map
>> >>> > to a
>> >>> > set of points which have non-integer coordinates and I want
>> the result
>> >>> > to be
>> >>> > as accurate as possible. Here is what I did, but the result
>> is not very
>> >>> > accurate.
>> >>> >
>> >>> > First I find the integer points which are the neighbors of
>> the target
>> >>> > points
>> >>> > and set these integer points as trial points. Then I use some
>> >>> > interpolation
>> >>> > method to initialize the distance from these trial points to
>> the target
>> >>> > points, which are assumed to be "exactly correct". The
>> TrialPoints in
>> >>> > the
>> >>> > FastMarchingImageFilter is defined as this set of trial
>> points and
>> >>> > their
>> >>> > corresponding distances to the target points. The AlivePoints
>> is empty.
>> >>> > When
>> >>> > I check the result distance map, I find that the distance
>> value of
>> >>> > these
>> >>> > trial points are changed, they are no longer what their
>> initial states
>> >>> > are.
>> >>> > Therefore, the iso curve deviate the original input a little
>> bit. I am
>> >>> > quite
>> >>> > confusing about this result.
>> >>> >
>> >>> > I noticed you mentioned on the mailing list about neighbor
>> update, that
>> >>> > is
>> >>> > to distinguish between user-specified trial points and
>> >>> > algorithm-generated
>> >>> > trial points. here is the discussion,
>> >>> > http://www.itk.org/pipermail/insight-users/2009-May/
>> 030282.html . I
>> >>> > wonder
>> >>> > if you have any suggestions about my problem.
>> >>> >
>> >>> > Any input is appreciated.
>> >>> >
>> >>> > Thanks
>> >>> > Siqi Chen
>> >>> >
>> >>
>> >
>> >
>> > _____________________________________
>> > 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
>> >
>> >
>>
>> _____________________________________
>> 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/20100108/06368a12/attachment-0001.htm>
More information about the Insight-users
mailing list