[Insight-users] Behavior of FastMarchingImageFilter

[Dan Mueller]

Hi Luca,

Sure -- adding the target stuff to the base class at the same time sounds great.

Do you want me to do it? I probably won't make the Friday 15th May
deadline for 3.14 though...

Cheers, Dan

2009/5/11 Luca Antiga <luca.antiga at gmail.com>:
> Hi Dan,
>  yep, that would work.
> The user-supplied Trial points should be tagged as InitialTrial in the
> LabelImage,
> put in the TrialHeap but not be updated as neighbors. And taking care of
> condition
>     if ( m_LabelImage->GetPixel( node.GetIndex() ) != TrialPoint )
> at line 292 in itkFastMarchingImageFilter.txx.
> Actually, since we are at it, it may also be good to move the Target
> business from
> FastMarchingUpwindGradientImageFilter to the base class (see your Apr 14
> message "Trying to monitor progression of fast march operation").
> Actually, we may also improve it by not looping around the target point list
> in
> the UpdateNeighbors method, but defining an extra label in the enum
> (Target)
> and storing the Target label in the LabelImage.
> I didn't do it back then because I just wanted to write an extra class
> without
> touching the original filter.
> What do you think?
>
> Luca
>
>
> On May 9, 2009, at 8:19 PM, Dan Mueller wrote:
>
> shouldn't you check if the update value computed for the trial point
>
> is smaller than the already existing value for that trial point and if
>
> so update it regardless?
>
> Ah yes, I see the issue with the fix I proposed -- once a trial point
> is given a solution value, it will never be improved upon.
>
> But what if the trial point is user specified? Shouldn't the value
> given by the user be used, even it is larger than the calculated one?
> Otherwise, why does the user even specify the node value?
>
> So I guess we need some way to distinguish between user specified and
> computed trial points... If the neighbour is a user specified trial
> point or alive, then don't update it...?
>
> Regards, Dan
>
> 2009/5/9 Luca Antiga <luca.antiga at gmail.com>:
>
> Hi Dan,
>
>  shouldn't you check if the update value computed for the trial point
>
> is smaller than the already existing value for that trial point and if
>
> so update it regardless?
>
> Luca
>
> 2009/5/9, Dan Mueller <dan.muel at gmail.com>:
>
> Hi Amar,
>
> Congratulations! It seems you managed to find a bug in the
>
> FastMarchingImageFilter :P
>
> In itkFastMarchingImageFilter.txx the UpdateNeighbors method checks if
>
> the candidate neighbour is an alive point. If so, then it is skipped
>
> from further processing. If not, then it is processed (ie. output
>
> value updated). This method should also test if the candidate
>
> neighbour is a *trial* point. If it is a trial point, then the output
>
> value has already been set, and does not need to be updated. Doing so
>
> will actually produce incorrect values, as you have encountered with
>
> your "jagged" edges.
>
> I have submitted a bug for this issue:
>
>     http://public.kitware.com/Bug/view.php?id=8990
>
> The bug report includes a patch, which I have also attached to this mail.
>
> Surprisingly this issue has been around since version 1.16 of the file
>
> (nearly 8 years)...
>
> Please let us know if this solves your issue.
>
> Regards, Dan
>
>
> 2009/5/8 Dan Mueller <dan.muel at gmail.com>:
>
> Hi Amar,
>
> I can reproduce your results. It is very strange; I can't see what is
>
> wrong with the code. I will look into it further on the weekend.
>
> Cheers, Dan
>
> 2009/5/8 Amardeep Singh <amar.singh at gmx.de>:
>
> Hello everybody:
>
> Dear Dan:
>
> Thank you very much for your help!
>
> I have attached some sample code and a binary image of a circle for
>
> testing.
>
> The program can be called like the following:
>
> fast_marching_itk_mailing_list ./circle.nii.gz
>
> ./fast_marching_output.nii.gz
>
> 200
>
> The last number is the stopping value which I usually set to large
>
> values. I
>
> use
>
> the fast marching filter in this example with a constant speed of 1 so
>
> that
>
> the result should be some kind of a distance map.
>
> Certainly, one important point is my calculation of the alive and trial
>
> points. All the points of the circle that have a background voxel within
>
> their
>
> 6-neighborhood are considered trial points (so the last layer around the
>
> circle becomes a layer of trial points). All other points of the circle
>
> are
>
> alive points.
>
> I have attached a screenshot of the 0 isoline of the result (the complete
>
> output image was too large for the mailing list, unfortunately). When you
>
> look at it, you can see the jagged edge again. But still: Shouldn't it
>
> propagate evenly from the trial points?
>
> The intialization (0 isoline of the result) looks jagged already.
>
> Concerning my previous snapshot:
>
> The output image was of type float. Just for the snapshopt, I mapped all
>
> values > 0 to 1 in order to visualize the jagged edge.
>
> Best regards
>
> Amar
>
>
> Dan Mueller wrote:
>
> Hi Amar,
>
> Perhaps you can post a minimal example (code, input image, output
>
> image) demonstrating the issue.
>
> I have not experienced this issue myself, but the
>
> FastMarchingImageFilter does propagate the front using
>
> 4-/6-connectivity, so after a small timestep it may be possible to
>
> experience the "jagged" edge you are seeing (though not the checked
>
> pattern inside the hole). How many iterations are you letting the
>
> filter run for? (ie. what is the stopping value?)
>
> Also, the FastMarchingImageFilter requires a real (float, double)
>
> output image. What is the pixel type of your output image?
>
> Cheers, Dan
>
> 2009/5/7 Amardeep Singh <amar.singh at gmx.de>:
>
>
> Hello everyone
>
> I am wondering about the behavior of the FastMarchingImageFilter.
>
> I want to let the algorithm march from a certain surface which is
>
> represented by a binary image. I
>
> pass the points of the object as alive points to the marcher. The layer
>
> outside of the object (in my case: a brain)
>
> is passed as trial points. When I look at the output of the fast
>
> marching filter, there is something strange -> see attached screenshot.
>
> The area within the yellow line is my initial object. So they are all
>
> alive points whose value is 0 in the output.  But when you look at the
>
> layer
>
> outside of the inner object, you see that it is jagged, black (value 0)
>
> and white (value != 0) voxels take turns. Why is this the case?
>
> What do I need to do, so that the marcher marches equally away from the
>
> alive points. The speed image values on the layer around the alive
>
> points do have the same value, so this is not the reason.
>
> Thanks a lot in advance!
>
> Best Regards
>
> Amar