[Insight-developers] level sets: curvature term and time step
Joshua Cates
cates at sci . utah . edu
Thu, 10 Jul 2003 16:03:04 -0600 (MDT)
Hi Paul and Lydia,
Yes I agree that the CurvatureSpeed() should be taken account in the
time step calculation. I'll have a look at this and get it resolved.
Josh.
______________________________
Josh Cates
School of Computer Science
University of Utah
Email: cates at sci . utah . edu
Phone: (801) 587-7697
URL: http://www . sci . utah . edu/~cates
On Thu, 10 Jul 2003, Paul Yushkevich wrote:
> I agree, in SNAP we use the curvature speed, so it would be better if it
> were taken into consideration.
>
> Paul
>
> Lydia Ng wrote:
>
> >Hi Josh,
> >
> >I think it is possible (but may not be likely) that all the level set
> >filters (ie. using itk::LevelSetFunction) are affected not just the
> >spatially varying curvature ones.
> >
> >We have,
> >Curvature_term = curvature * CurvatureWeight * CurvatureSpeed();
> >
> >
> >Say for example, there is no advection term then in
> >ComputeGlobalTimeStep
> >
> >dt = m_DT / vnl_math_abs( m_CurvatureWeight );
> >
> >
> >By default, CurvatureSpeed() returns 1 - so things might go wrong if
> >curvature is much larger than 1 / m_DT ( > 4 for 2D ).
> >
> >But I don't know how likely this scenario is.
> >
> >________________
> >
> >The problem Nils ( email attached ) is having is because his
> >CurvatureSpeed() returns 1000 and he has a CurvatureWeight of 1.
> >
> >Say if curvature is approx 1 then curvature_term is 1000
> >But dt = m_DT = 0.25 (2D) and hence the algorithm takes to big a step.
> >
> >While if he had it the other way around: CurvatureSpeed() = 1 and
> >CurvatureWeight = 1000 then dt = 0.25/1000 and the algorithm remains
> >stable.
> >
> >So I think the time step calculation to be general it needs to take into
> >account the whole curvature_term rather than just the CurvatureWeight.
> >
> >- Lydia
> >
> >
> >
> >>-----Original Message-----
> >>From: Joshua Cates [mailto:cates at sci . utah . edu]
> >>Sent: Thursday, July 10, 2003 10:32 AM
> >>To: Lydia Ng
> >>Cc: Insight-Developers (E-mail)
> >>Subject: Re: [Insight-developers] level sets: curvature term and time
> >>
> >>
> >step
> >
> >
> >>Hi Lydia,
> >>
> >>Just to clarify: you are talking about the spatially varying curvature
> >>term, and not the curvature calculation itself? So only filters
> >>which use this term should be affected, right?
> >>
> >>Josh.
> >>
> >>______________________________
> >> Josh Cates
> >> School of Computer Science
> >> University of Utah
> >> Email: cates at sci . utah . edu
> >> Phone: (801) 587-7697
> >> URL: http://www . sci . utah . edu/~cates
> >>
> >>
> >>On Thu, 10 Jul 2003, Lydia Ng wrote:
> >>
> >>
> >>
> >>>Question for those interested in level sets:
> >>>
> >>>Currently the itk::LevelSetFunction class does not take the
> >>>
> >>>
> >curvature
> >
> >
> >>>term into account when computing the time step to meet the CFL
> >>>
> >>>
> >criterion
> >
> >
> >>>to ensure numerical stability.
> >>>
> >>>This is causing some grief with users playing around with the
> >>>
> >>>
> >scaling
> >
> >
> >>>parameter in each of the term in particular when they set the
> >>>
> >>>
> >curvature
> >
> >
> >>>scaling to be high relative to the propagation and advection terms.
> >>>
> >>>Can someone point me to a paper/book that talks about have to
> >>>
> >>>
> >compute
> >
> >
> >>>the CFL time step for the curvature term?
> >>>
> >>>Alternatively, (hacking versus mathematically sound) could we just
> >>>restrict the time step such that the change in level set value due
> >>>
> >>>
> >to
> >
> >
> >>>the curvature term is no larger that half a pixel? Would this have
> >>>
> >>>
> >the
> >
> >
> >>>desired affect?
> >>>
> >>>The second alternative is for me to change my filters
> >>>ShapeDetectionLevelSetImageFilter and
> >>>GeodesicActiveContourLevelSetImageFilter just to use a constant time
> >>>step - making it the user's problem to ensure they choose the time
> >>>
> >>>
> >step
> >
> >
> >>>small enough. This is how it was before, when I got complaints from
> >>>other users about how it was difficult to set this parameter...
> >>>
> >>>Cheers,
> >>>Lydia
> >>>_______________________________________________
> >>>Insight-developers mailing list
> >>>Insight-developers at itk . org
> >>>http://www . itk . org/mailman/listinfo/insight-developers
> >>>
> >>>
> >>>
> >
> >
> >
> >
> > ------------------------------------------------------------------------
> >
> > Subject:
> > RE: [Insight-users] Crash in ShapeDetectionLevelSetImageFilter
> > From:
> > "Nils Hanssen" <hanssen at caesar . de>
> > Date:
> > Thu, 10 Jul 2003 00:16:22 -0700
> > To:
> > "Luis Ibanez" <luis . ibanez at kitware . com>
> >
> >
> >Hi Luis,
> >
> >thank you for the detailed answer.
> >
> >
> >My initial level set conforms to the "negative inside" convention.
> >Now, I am using a curvature weighting of 1000 (just to see what happens in
> >the extreme case) and a propagation weighting of 0.001. When running the
> >ShapeDetectionLSIF, the contour is not collapsing at all. Maybe that's
> >because of the feature image that dictates the curvature speed (==
> >propagation speed in the ShapeDetectionLSIF)?
> >
> >Since the curvature speed in the ShapeDetectionLSF equals the propagation
> >speed, I derived my own classes called "TubularLevelSetFunction" and
> >"TubularLevelSetImageFilter". These are basically clones of the
> >ShapeDetection* classes, but with the difference that I set the curvature
> >speed to (constant) one in the Level-set function.
> >Now, the contour is collapsing according to the mean curvature but all this
> >happens _very_ slowly (the curvature weighting is 1000!).
> >
> >In itk::LevelSetFunction::ComputeUpdate(...) I see the following:
> >---
> > curvature_term = curve;
> > curvature_term *= m_CurvatureWeight * this->CurvatureSpeed(it, offset);
> >---
> >so it should make no difference when I set the curvature weight to 1 and the
> >curvature speed to (constant) 1000 (and not vice versa), right? But when I
> >do this, I get numerical instabilities during the evolution. What is wrong
> >with my assumption?
> >
> >Luis, what was the curvature weighting and curvature speed when you tested
> >the collapsing contour with the syntethic sphere image? Was the curvature
> >speed constant or dependent on the image?
> >
> >For now, I do all calculations in 2D. Could that be a problem? I have
> >installed a observer that shows me the progress of the contour after each
> >update event.
> >
> >
> >Regards,
> >Nils
> >
> >
> >
> >
> >>-----Original Message-----
> >>From: Luis Ibanez [mailto:luis . ibanez at kitware . com]
> >>Sent: Wednesday, July 09, 2003 6:14 PM
> >>To: Nils Hanssen
> >>Cc: Insight-users at public . kitware . com
> >>Subject: Re: [Insight-users] Crash in
> >>ShapeDetectionLevelSetImageFilter
> >>
> >>
> >>
> >>Hi Nils,
> >>
> >>You are right in your interpretation of the
> >>behavior for the ShapeDetectionLevelSetFilter.
> >>
> >>Increasing the curvature scaling while decreasing
> >>the propagation scaling should make the contour
> >>collapse (contract).
> >>
> >>I just verified this behavior using a syntethic
> >>image of a sphere.
> >>
> >>Note that in version 1.2 this filter is expecting
> >>the input level set to conform to the convention
> >>of "negative inside" which means that the level
> >>set has negative values 'inside' the contour and
> >>positive values outside the contour.
> >>
> >>You may have to push the curvature scaling to
> >>a value higher than 1.0. (e.g. 5.0 or so).
> >>
> >>I would avoid to use 0.0 for the propagation
> >>scaling ( in part just for superstition against
> >>zeros values as parameters). You may want to
> >>try something like 0.1 for the propagation
> >>scaling.
> >>
> >>Please verify the convention used by your initial
> >>level set. That may be the cause for the contour
> >>not behaving as expected.
> >>
> >>
> >>Regards,
> >>
> >>
> >> Luis
> >>
> >>
> >>
> >>---------------------
> >>Nils Hanssen wrote:
> >>
> >>
> >>>Hi Luis,
> >>>
> >>>I am using the 1.2.0 web release.
> >>>
> >>>Regards,
> >>>Nils
> >>>
> >>>
> >>>
> >>>
> >>>>-----Original Message-----
> >>>>From: insight-users-admin at itk . org
> >>>>[mailto:insight-users-admin at itk . org]On
> >>>>Behalf Of Luis Ibanez
> >>>>Sent: Tuesday, July 08, 2003 3:55 PM
> >>>>To: Nils Hanssen
> >>>>Cc: Insight-users at public . kitware . com
> >>>>Subject: Re: [Insight-users] Crash in
> >>>>ShapeDetectionLevelSetImageFilter
> >>>>
> >>>>
> >>>>
> >>>>Hi Nils,
> >>>>
> >>>>Are you using a very recent checkout of ITK ?
> >>>>
> >>>>The level set filters have been reworked in
> >>>>recent days, so it will help us to know if you
> >>>>are experiencing this behavior with the current
> >>>>version or with the version as it was two weeks
> >>>>ago (or before that).
> >>>>
> >>>>Please let us know,
> >>>>
> >>>> Thanks
> >>>>
> >>>>
> >>>> Luis
> >>>>
> >>>>
> >>>>------------------------
> >>>>Nils Hanssen wrote:
> >>>>
> >>>>
> >>>>
> >>>>>Hi,
> >>>>>
> >>>>>I am trying to understand the behaviour of the ShapeDetectionLSIF.
> >>>>>Therefore, I set the propagation-weighting to zero. By setting the
> >>>>>curvature-weighting to a value of one I would expect that
> >>>>>
> >>>>>
> >>>>the inital
> >>>>
> >>>>
> >>>>
> >>>>>surface is shrinking to a point (I set the MaxRMSError to
> >>>>>
> >>>>>
> >>>>zero) and the
> >>>>
> >>>>
> >>>>
> >>>>>number of iterations very high.
> >>>>>
> >>>>>The filter is crashing in
> >>>>>SegmentationLevelSetFunction<TImageType, TFeatureImageType>
> >>>>>::PropagationSpeed(const NeighborhoodType &neighborhood, const
> >>>>>FloatOffsetType &offset) const
> >>>>>[...]
> >>>>>--> else return (
> >>>>>
> >>>>>
> >>>>static_cast<ScalarValueType>(m_SpeedImage->GetPixel(idx)) );
> >>>>// crashing here
> >>>>
> >>>>
> >>>>
> >>>>>[...]
> >>>>>
> >>>>>When I set the propagation-weighting to - for instance -
> >>>>>
> >>>>>
> >>0.0001 the
> >>
> >>
> >>>>>filter is not crashing, but the contour is not shrinking
> >>>>>
> >>>>>
> >>>>according to
> >>>>
> >>>>
> >>>>
> >>>>>the mean curvature.
> >>>>>
> >>>>>Is that the correct behaviour of the filter?
> >>>>>
> >>>>>I would appreciate any suggestions. Thanks!
> >>>>>
> >>>>>
> >>>>>Regards,
> >>>>>Nils
> >>>>>
> >>>>>
> >>>>>-------------------------
> >>>>>Nils Hanssen
> >>>>>Surgical Systems Laboratory
> >>>>>research center c ae sa r
> >>>>>Ludwig-Erhard-Allee 2
> >>>>>53175 Bonn
> >>>>>fon: +49-228-9656-197
> >>>>>fax: +49-228-9656-117
> >>>>>___http://www . caesar . de/ssl_
> >>>>>
> >>>>>
> >>>>>
> >>>>>
> >>>>>
> >>>>
> >>>>_______________________________________________
> >>>>Insight-users mailing list
> >>>>Insight-users at itk . org
> >>>>http://www . itk . org/mailman/listinfo/insight-users
> >>>>
> >>>>
> >>>>
> >>>_______________________________________________
> >>>Insight-users mailing list
> >>>Insight-users at itk . org
> >>>http://www . itk . org/mailman/listinfo/insight-users
> >>>
> >>>
> >>>
> >>
> >>
> >>
> >>
> >
> >_______________________________________________
> >Insight-users mailing list
> >Insight-users at itk . org
> >http://www . itk . org/mailman/listinfo/insight-users
> >
> >
> >
> >
>
>
> --
> --------------------------------
> Paul A. Yushkevich, Ph.D.
> President, Cognitica Corporation
>
> 17 Flemington Rd
> Chapel Hill, NC 27517
> Tel: 1-919-929-7652
> --------------------------------
>
>
>