[Insight-users] Behavior of FastMarchingImageFilter

[Dan Mueller]

Hi Amardeep,

The initial patch I proposed missed a use case, so: no, the results
you experienced are not normal. The solution Luca proposed (adding an
InitialTrial type) should correctly solve your issue.

I'm snowed under at work at the moment, so I won't get around to
submitting the fix until the weekend...

Regards, Dan

2009/5/12 Amardeep Singh <amar.singh at gmx.de>:
> Hello
>
> I would just like to let you know what I have experienced with the bugfix as
> it is
> now.
> The "jaggednes" of the result has definitely disappeared. However, I am
> confused
> by the fact that the isolines of the fast marching output a little further
> away from the object look like edges of a square
> rotated by 45 degrees (see the snapshot at the isoline 12 and the attached
> output image; I have applied the filter to the
> binary circle again).
> Is this normal? I understand that the current bugfix is a preliminary one.
>
> I am looking forward to the final bugfix and thanks a lot for correcting
> this bug.
>
> Best regards
> Amar
>
>
> Dan Mueller wrote:
>>
>> 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
>>>
>>
>>
>