[Insight-users] Concatenation of deformation fields?
Tom Vercauteren
tom.vercauteren at m4x.org
Tue Apr 6 10:05:01 EDT 2010
Hi Gerald,
On Tue, Apr 6, 2010 at 15:23, Lodron, Gerald <Gerald.Lodron at joanneum.at> wrote:
> I read your paper and i find it really interresting since it calculates the inverse and the forward transformation/deformation and have a few questions, may you can help me understanding:
It's nice to hear that my work may help you.
> - In the paper you have the terms "spatial transformation (a vector field) s" and the term "velocity field" which also seems to be a vector field. What is the difference (interpretation would be interesting, the formulation with exp() seems a little bit strange for me, sorry for my bad mathematics and image processing skills) between them?
I'll try to clarify things a little bit. If you need more than this
hand-waving explanations, please refer to:
Arsigny et al., MICCAI 2006
http://www.inria.fr/sophia/asclepios/Publications/Arsigny/arsigny_miccai06.pdf
and
Vercauteren et al., NeuroImage 2009
http://www.inria.fr/sophia/asclepios/Publications/Tom.Vercauteren/DiffeoDemons-NeuroImage08-Vercauteren.pdf
and
Vercauteren et al., MICCAI 2008
http://www.inria.fr/sophia/asclepios/Publications/Tom.Vercauteren/SymLogDemons-MICCAI08-Vercauteren.pdf
Let p be a point and v a deformation field. If you want to map p
through v you will end up with a point q whose coordinates are:
q = p + v(p)
Indeed the deformation field only provides deformation with respect to
the identity transformation (i.e. no deformation maps a point to
itself). The spatial transformation T corresponding to the deformation
field v is thus given by
T = Id + v
As illustrated by the equations, a deformation field provides for each
point the complete deformation vector for the transformation. Hence it
is represented by a vector field.
A velocity field is something a little more subtle. You have to think
in terms of kinematics. Instead of saying for each point where it will
end up, you only provide the direction of its motion and its speed,
i.e. you provide a velocity vector for each point. A velocity field is
a vector field that gives this velocity information for each point in
space. If you want to compute the deformation from the velocity you
will have to integrate the velocity over time. The exponential I am
referring to makes the simple assumption that the velocity field does
not change over the course of time and integrates it over a unit
amount of time. With this exponential notation, the spatial
transformation S corresponding to a velocity field w is given by
S = exp(w)
If you need the deformation field that corresponds to this
transformation, you will need to subtract the identity, i.e. if u is
the deformation field corresponding to S, you will have
S = exp(w) = Id + u
u = exp(w) - Id
or for each point
u(p) = exp(w)(p) - p
Note that it is always possible to represent an exponential field with
a deformation field but the converse is not always true.
> - What is the difference between itk::TransformToDeformationFieldSource and itk::TransformToVelocityFieldSource ?
TransformToDeformationFieldSource gives v when you have T with a parametric form
TransformToVelocityFieldSource gives w when you have S with a parametric form
Hope this helps,
Tom
> -----Ursprüngliche Nachricht-----
> Von: Tom Vercauteren [mailto:tom.vercauteren at gmail.com]
> Gesendet: Samstag, 27. März 2010 19:30
> An: Lodron, Gerald
> Cc: Luis Ibanez; itk
> Betreff: Re: [Insight-users] Concatenation of deformation fields?
>
> Hi Gerald,
>
> This process might actually be eased by a utility class called
> DisplacementFieldCompositionFilter
> that you will find in the following IJ paper:
> http://www.insight-journal.org/browse/publication/644
>
> This filter does exactly the combination of steps 2 and 3 proposed by Luis.
>
> Hope this helps,
> Tom
>
> On Sat, Mar 27, 2010 at 17:04, Luis Ibanez <luis.ibanez at kitware.com> wrote:
>> Hi Gerald,
>>
>> Thanks for your very detailed clarification.
>>
>> Just to simplify the email notation,
>> let me add the following two definitions:
>>
>> A := D_r^-1
>> B := D_n^-1
>> C := D_n,r^-1
>>
>>
>> so the deformation field "C" is the result of composing the
>> deformation field "A"
>> with "B".
>>
>> Therefore, for every point x, its deformed target x" will be computed
>> as:
>>
>> x" = C(x) = B( A (x) )
>>
>> that is
>>
>> x' = A( x )
>> x" = B( x' )
>>
>> x" = C( x ) = B( x' ) = B( A ( x ) )
>>
>>
>> Let's now name the coordinate frames:
>>
>> F1 := The original image from
>> which the points "x" are taken.
>> This is the undeformed image.
>>
>> F2 := The deformed image that results
>> from applying the field A to all the
>> points of F1.
>>
>> F3 := The deformed image that results
>> from applying the field B to all the
>> points of F2.
>>
>>
>> With your current processing you have
>> direct access to the fields "A" and "C".
>>
>> "A" maps the coordinate system F1 to F2
>> "C" maps the coordinate system F1 to F3
>>
>> and you are interested in finding the field B where
>>
>> "B" maps the coordinate system F2 to F3
>>
>> So one way to *estimate* the field "B" is the following.
>>
>> 1) Use the filter
>>
>> Insight/Code/BasicFilter/
>> itkIterativeInverseDeformationFieldImageFilter.h
>>
>> in order to estimate A^-1, (let's call it A' so it looks
>> better in an email: A' = A^-1.
>>
>> A' maps the coordinates from F2 to F1
>>
>> 2) Having A' you can now visit all the points x' of F2,
>> map them to x in F1 via A', and then map their
>> destinations to F3 via C.
>>
>> x'' = C ( A'( x' ) )
>>
>> Therefore you can build an estimation of B by
>> composing C with A'
>>
>> 3) Note that when you apply A' to a point x', the
>> resulting point x will not land exactly in the grid
>> node of the image F1, and therefore, you will
>> have to interpolate from C in order to find the
>> corresponding destination.
>>
>> For this purpose you can use the class
>>
>> Insight/Code/Common/
>> itkVectorLinearInterpolateImageFunction.h
>>
>>
>>
>> If you compose steps (1),(2),(3) in to an ITK filter, that will make a
>> very nice contribution to the Insight Journal :-)
>>
>>
>>
>> Regards,
>>
>>
>> Luis
>>
>>
>> ----------------------------------------------------------------------
>> ------------------------------ On Fri, Mar 19, 2010 at 3:30 AM,
>> Lodron, Gerald <Gerald.Lodron at joanneum.at> wrote:
>>> Hi
>>>
>>> I do not exactly know what you mean but i will try to explain my pipline in a little more detail:
>>>
>>> In principle i have two 3D data sets with I1 = 512 x 512 x z1 and I2 = 512 x 512 x z2. The origin and pixel spacing could be different. So the first thing i do is to calculate a new image I2iso = f(I1) and I2iso = f(I2) so that the pixel spacing becomes isotrophic, i mean pixel spacing in x y and z is now the same and also equal for I2iso and I1iso (spacingX1=spacingY1=spacingZ1=spacingX2=spacingY2=spacingZ2). I need an isotrophic pixel spacing because i use the demon registration in my nonrigid step (and the usepixelspacing function of the filter always leads to an error).
>>>
>>>
>>> Now i made two registrations, a coarse one (i called it rigid which is wrong) and a fine one where the coarse one is implemented via two subregistrations. First i make a versor registration, so i compute translation and rotation. Then i make a coarse Bspline like in the examples with 4 grid points and the versor result as initial parameters. The result is converted into a deformation field (for sure in reference to the fixed image), this is called my D_r^-1 (i know it is not pure rigid, but nearly).
>>>
>>> Now i make my nonrigid registration using a multiscale demon registration filter with the initial deformation field of the quasi rigid step (D_r^-1). So the result is the combination of the known D_r^-1 and the unknown D_n^-1 which is called D_n,r^-1.
>>>
>>> I want to display changes (which should refer to the norm of D_n-1) in the registered image (so I1 transformed with D_n,r^-1) so i need to extract D_n-1 from D_n,r^-1 using D_r^-1 and thats not a simply substraction...
>>>
>>>
>>>
>>>
>>>
>>> -----Ursprüngliche Nachricht-----
>>> Von: Luis Ibanez [mailto:luis.ibanez at kitware.com]
>>> Gesendet: Donnerstag, 18. März 2010 17:32
>>> An: Lodron, Gerald
>>> Cc: insight-users at itk.org
>>> Betreff: Re: [Insight-users] Concatenation of deformation fields?
>>>
>>> Hi Gerald,
>>>
>>> It should be possible to separate the non rigid component from the final deformation field.
>>>
>>> Did you compute both the
>>>
>>> 1) D_r^1 and
>>> 2) D_n^1
>>>
>>> in the same image grid ?
>>>
>>> Note however that the key point here is a problem of reference systems. The non-ridig deformation field will be one with respect to the rotated image, not with respect to the original image.
>>>
>>>
>>> Luis
>>>
>>>
>>> ---------------------------------------------------------------------
>>> ----------------- On Wed, Mar 10, 2010 at 11:52 AM, Lodron, Gerald
>>> <Gerald.Lodron at joanneum.at> wrote:
>>>>
>>>> Hi,
>>>>
>>>> I have a quite interresting problem and I am not sure if it is solveable but maybe someone can help:
>>>>
>>>> I made a rigid registration and converted the output to a deformation field, let us name it D_r^-1 (it is the inverse deformation because of inverse mapping). After that i made a nonrigid registration with the R_rigid^-1 as start deformation value, the output is the inverse deformation field of the rigid and the nonrigid registration, let us call it D_r,n^-1.
>>>>
>>>> My question is, is it possible to calculate the pure nonrigid deformation field D_n^-1 ?
>>>>
>>>> Best regards
>>>>
>>>>
>>>> _____________________________________
>>>> Powered by www.kitware.com
>>>>
>>>> Visit other Kitware open-source projects at
>>>> http://www.kitware.com/opensource/opensource.html
>>>>
>>>> Kitware offers ITK Training Courses, for more information visit:
>>>> http://www.kitware.com/products/protraining.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
>>
>> Kitware offers ITK Training Courses, for more information visit:
>> http://www.kitware.com/products/protraining.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
>>
>
More information about the Insight-users
mailing list