[Insight-users] Advection in ThresholdLevelSetImageFilter
Luca Antiga
luca.antiga at gmail.com
Tue Sep 23 10:04:30 EDT 2008
Dear Dário,
you're right, if you build a vector image and you advect it, yes,
it should
be treated as an advection field. Be careful though: if you compute a
vector,
like a gradient vector, out of a Laplacian (if I understood your
point correctly),
you're going to rely on 3rd derivatives: noise is going to kill you
at that point.
Even the laplacian itself is very sensitive to noise: if you don't 1.
initialize
your level sets close to the expected edges and/or 2. adopt a smooth
kernel
to compute the derivatives, you're not going to get good results.
With 3rd derivatives, you're going to make things even more severe,
so that
you may end up smoothing your data too much for your level set to be
usable.
I'm not saying it won't work, but just so you know.
As to curvature, the \kappa is actually the formulation of the mean
curvature
of a level set. \kappa is computed in the ComputeCurvatureTerm method in
itk::LevelSetFunction. The actual level set change due to curvature
is computed as
m_CurvatureWeight * this->CurvatureSpeed(...) * this-
>ComputeCurvatureTerm(...)
where the CurvatureSpeed is the Z you're looking for.
If you look at itkGeodesicActiveContourLevelSetFunction.h, you'll see
that CurvatureSpeed
is
/** The curvature speed is same as the propagation speed. */
virtual ScalarValueType CurvatureSpeed(const NeighborhoodType &
neighborhood,
const FloatOffsetType &
offset, GlobalDataStruct *gd ) const
{ return PropagationSpeed( neighborhood, offset, gd ); }
so, basically, the smaller the speed, the smaller is the effect of
curvature, so as to
say, you regularize where the feature content is low, while you
regularize less close
to you real features (remember that propagation speed is low at
features, opposite
to the gradient magnitude).
And this is actually the opposite to what you're suggesting :-), but
it gives you an example
on how to derive a custom Function that does it the way you want it.
I hope I addressed your questions, but let me know.
Cheers
Luca
On Sep 23, 2008, at 3:42 PM, Dário Oliveira wrote:
> Dear Luca, thank you.
>
> If I understood properly what you've said, I agree with you. Could I
> build a vector field using the laplacian 'image'? I mean, in this case
> the vector field would have high absolute divergence values near
> edges, and null values exactly at the edges, and then hopefully it
> would atract the interface to these points, and then lock there. If
> so, then this term would be 'usable' in Advection, as it was defined
> in level sets theory. Is that right?
>
> Let me ask you one more question? In the curvature term, I have this Z
> and a \kappa. This \kappa is computed as the divergence of the unit
> normal vector to the interface, where high values indicate high
> curvature, and low values low curvature. But this Z image, I don't
> figured out. I can set this image? It should be something like the
> gradient of my input image? In this case, I would only get non-null
> curvature values near edges, is that correct?
>
> Sorry for so many questions.
>
> Thanks,
> Dário
>
> 2008/9/23 Luca Antiga <luca.antiga at gmail.com>:
>> Dear Dario,
>> I must say I don't agree with you.
>> "Advection" means that the level set term is of type \grad P \cdot
>> \grad
>> \phi,
>> where P is an advection potential, such as | \grad I |, I being
>> the image,
>> as in
>> geodesic active contours.
>> The way Laplacian zero-crossing acts is instead through the
>> "propagation"
>> term, G | \grad \phi |, where G is a proper advection speed, which
>> in case
>> of
>> Laplacian zero crossings is simply the value of the laplacian
>> (which is zero
>> at the gradient magnitude ridges).
>>
>> So, from a "purpose" point of view you're right, the laplacian
>> tends to
>> attract
>> the surface to the ridges of the gradient magnitude, like the
>> advection term
>> does. However, it doesn't do it by advecting the level set with a
>> vector
>> field,
>> but by propagating the level set using the laplacian. That's why the
>> Laplacian
>> should stay where it is and not be moved to
>> CalculateAdvectionImage as
>> you're suggesting.
>>
>> Best regards
>>
>>
>> Luca
>>
>>
>>
>> On Sep 22, 2008, at 7:43 PM, Dário Oliveira wrote:
>>
>>> Thank you Luis!
>>>
>>> I have seen it yesterday. I notice that I can also reimplement the
>>> method CalculateAdvectionImage and then the advection term makes
>>> sense. I got confused because the way I understood the method, this
>>> part of calculating the zero-crossing of the Laplacian should be
>>> exactly in the CalculateAdvectionImage method. Actually, as far as I
>>> understood, it modifies the speed image in such a way that the
>>> interface locks at edges. I guess I will modify it to use the
>>> laplacian zero-crossings in the advection term.
>>>
>>> Thank you.
>>> Dário Oliveira
>>>
>>> 2008/9/22 Luis Ibanez <luis.ibanez at kitware.com>:
>>>>
>>>> Hi Dario,
>>>>
>>>> Yes, the EdgeWeight plays in the ThresholdSegmentationLevelSet
>>>> a simlar role to the AdvectionWeight in the GeodesicActiveContour.
>>>>
>>>> When you set the EdgeWeight to non-zero, the filter computes a
>>>> Laplacian and uses it to drive the zero set closer to the
>>>> zero-crossings of the Laplacian.
>>>>
>>>>
>>>> The SmoothingConductance is used as a parameter for the Gradient
>>>> Anisotropic diffusion filter that is used inside of the Threshold
>>>> SegmentationLevelSetFunction, to smooth the image before computing
>>>> the Laplacian on it. Otherwise, the Laplacian will find a large
>>>> number of zero-crossings in all the small details (high spatial
>>>> frequency) of the image.
>>>>
>>>>
>>>> BTW: You can review all this by simply looking at the source code
>>>> of the file:
>>>>
>>>>
>>>> Insight/Code/Algorithms/
>>>> itkThresholdSegmentationLevelSetFunction.txx
>>>>
>>>>
>>>> in particular by looking at
>>>>
>>>>
>>>> lines 43-57:
>>>> =============
>>>>
>>>> if (m_EdgeWeight != 0.0)
>>>> {
>>>> diffusion->SetInput(this->GetFeatureImage());
>>>> diffusion->SetConductanceParameter(m_SmoothingConductance);
>>>> diffusion->SetTimeStep(m_SmoothingTimeStep);
>>>> diffusion->SetNumberOfIterations(m_SmoothingIterations);
>>>>
>>>> laplacian->SetInput(diffusion->GetOutput());
>>>> laplacian->Update();
>>>>
>>>> lit = ImageRegionIterator<FeatureImageType>(laplacian-
>>>> >GetOutput(),
>>>>
>>>> this->GetFeatureImage()->GetRequestedRegion());
>>>> lit.GoToBegin();
>>>> }
>>>>
>>>>
>>>> and lines 77-85:
>>>> ================
>>>>
>>>>
>>>> if ( m_EdgeWeight != 0.0)
>>>> {
>>>> sit.Set( static_cast<ScalarValueType>(threshold +
>>>> m_EdgeWeight *
>>>> lit.Get()) );
>>>> ++lit;
>>>> }
>>>> else
>>>> {
>>>> sit.Set( static_cast<ScalarValueType>(threshold) );
>>>> }
>>>>
>>>>
>>>> ----------------
>>>>
>>>>
>>>> "Glimpsing at the Source leaves no Doubt"
>>>>
>>>>
>>>>
>>>> Regards,
>>>>
>>>>
>>>> Luis
>>>>
>>>>
>>>> -----------------------
>>>> Dário Oliveira wrote:
>>>>>
>>>>> Dear all,
>>>>>
>>>>> I guess I miss something here. Why the advection scaling term
>>>>> is not
>>>>> applied in ThresholdLevelSetImageFilter? I performed some
>>>>> tests, and
>>>>> the value of advection scaling term makes no difference at all
>>>>> in the
>>>>> final result. Does the EdgeWeight term play the role of advection
>>>>> scaling?
>>>>>
>>>>> Another question... What does the smoothing conductance term?
>>>>>
>>>>> Thanks for your attention!
>>>>
>>>
>>>
>>>
>>> --
>>> Dário Oliveira
>>> _______________________________________________
>>> Insight-users mailing list
>>> Insight-users at itk.org
>>> http://www.itk.org/mailman/listinfo/insight-users
>>
>>
>
>
>
> --
> Dário Oliveira
--
Luca Antiga, PhD
Head, Medical Imaging Unit,
Bioengineering Department,
Mario Negri Institute
email: antiga at marionegri.it
web: http://villacamozzi.marionegri.it/~luca
mail: Villa Camozzi, 24020, Ranica (BG), Italy
phone: +39 035 4535-381
More information about the Insight-users
mailing list