[Insight-users] Computing the gradient in the Metric?

Luis Ibanez luis.ibanez at kitware.com
Wed Sep 9 12:18:30 EDT 2009


Hi Motes


Let's call

x = physical point in the fixed image coordinate system

T = transform

p = Transform parameters

x' = T(p,x) = physical point x mapped to moving image
       coordinate system by using the Transform T with
       parameters p.

J = dx' / dp  = Jacobian of mapped point x' with respect
       to the Transform parameters p.


Therefore, what Luke meant is that
you can use the chain rule as:

   dM / dp = ( dM / dx' ) * ( dx' / dp )

where

dM / dx'  = Gradient vector of the Moving image

and

dx' / dp = Jacobian of the Transform.



     Regards,


           Luis


------------------------------------------------------------------------------
On Sun, Sep 6, 2009 at 3:38 PM, motes motes <mort.motes at gmail.com> wrote:

> Hm not sure I understand:
>
>    dM/dx(T(p,x))  * dT/dp:
>
>
>
> Don't you mean:
>
>    dM/dx'  * dT(p,x)/dp:
>
> where:
>
>    x' = T(p,x)
>
> corresponding to the transformed pixel in the moving image M?
>
>
>
>
>
> On Sun, Sep 6, 2009 at 5:07 PM, Luke Bloy<luke.bloy at gmail.com> wrote:
> > You use the chain rule.
> >
> > dt/dp = dM(T(p,x))/dp = dM/dx(T(p,x)) * dT/dp
> >
> > so dM/dx(T(p,x)) is the gradient of the moving image evaluated at T(p,x)
> and
> > dT/dp is the jacobian of the transform.
> >
> > -luke
> >
> > motes motes wrote:
> >>
> >> Ok I am trying to do the metric E differentiation by hand:
> >>
> >> E = \sum[F(x) - M(T(p,x))]^2
> >>   = t^2
> >>
> >> where:
> >>
> >> t = F(x) - M(T(p,x))
> >>
> >> dE/dp = 2t * dt/dp
> >>
> >> where:
> >>
> >> dt/dp = dM(T(p,x))/dp
> >>
> >> since the fixed image is treated as a constant because we are
> >> differentiating with respect to the parameters p.
> >>
> >>
> >> Now I think I understand why its the gradient of the moving image that
> >> is used. But how are the Jacobian introduced in:
> >>
> >> dt/dp = dM(T(p,x))/dp
> >>
> >> ?
> >>
> >>
> >>
> >>
> >> On Sun, Sep 6, 2009 at 3:30 PM, motes motes<mort.motes at gmail.com>
> wrote:
> >>
> >>>
> >>> Thanks for the hint I have read those pages but I still have two
> >>> questions to the below expression:
> >>>
> >>>         sum += 2.0 * diff * jacobian( dim, par ) * gradient[dim];
> >>>
> >>> 1) In the above expression gradient[dim] corresponds to the gradient
> >>> for the current pixel in the moving image? On page 249 in the "Insight
> >>> Into Images" this is the first term in (10.6). But why are the
> >>> gradient of the moving image used to compute the derivative of the
> >>> metric? I know this is a very basic question but I still cannot see
> >>> the connection.
> >>>
> >>> 2) When I debug the code the jacobian actually just contains the
> >>> weights (computed using the transform kernel) for each node in the
> >>> support region. But are the jacobian not supposed to contain the first
> >>> order partial derivatives? It seems a bit confusing that the jacobian
> >>> is used as a container for the kernel transform weights.
> >>>
> >>>
> >>>
> >>>
> >>>
> >>> On Sun, Sep 6, 2009 at 1:01 PM, Neuner Markus<neuner.markus at gmx.net>
> >>> wrote:
> >>>
> >>>>
> >>>> Hi,
> >>>>
> >>>> To understand why the jacobian is used i would suggest to read pages
> >>>> 249-251
> >>>> in the Book "Insight Into Images" by Terry S. Yoo. Ther eis described
> >>>> how
> >>>> and why the jacobian is used.
> >>>>
> >>>> Greets
> >>>>
> >>>> motes motes wrote:
> >>>>
> >>>>>
> >>>>> In itkMeanSquaresImageToImageMetric.txx the metric value and the
> >>>>> gradient is computed. I pretty much understand how the metric value
> is
> >>>>> computed but am a bif confused on how the gradient is computed:
> >>>>>
> >>>>>
> >>>>>
> >>>>>
> >>>>>     const RealType movingValue  = this->m_Interpolator->Evaluate(
> >>>>> transformedPoint );
> >>>>>     const TransformJacobianType & jacobian =
> >>>>> this->m_Transform->GetJacobian( inputPoint );
> >>>>>     const RealType fixedValue     = ti.Value();
> >>>>>     this->m_NumberOfPixelsCounted++;
> >>>>>     const RealType diff = movingValue - fixedValue;
> >>>>>     measure += diff * diff;
> >>>>>
> >>>>>     for(unsigned int par=0; par<ParametersDimension; par++)
> >>>>>       {
> >>>>>       RealType sum = NumericTraits< RealType >::Zero;
> >>>>>       for(unsigned int dim=0; dim<ImageDimension; dim++)
> >>>>>         {
> >>>>>         sum += 2.0 * diff * jacobian( dim, par ) * gradient[dim];
> >>>>>         }
> >>>>>       derivative[par] += sum;
> >>>>>       }
> >>>>>     }
> >>>>>
> >>>>> It basically comes down to this line:
> >>>>>
> >>>>>
> >>>>>
> >>>>>         sum += 2.0 * diff * jacobian( dim, par ) * gradient[dim];
> >>>>>
> >>>>> why is the jacobian multiplied with the current gradient?
> >>>>> _____________________________________
> >>>>> Powered by www.kitware.com
> >>>>>
> >>>>> Visit other Kitware open-source projects at
> >>>>> http://www.kitware.com/opensource/opensource.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
> >>
> >> 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
>
> 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/20090909/c3035226/attachment.htm>


More information about the Insight-users mailing list