[slicer-devl]Tensor scalar calculation in Python

classic Classic list List threaded Threaded
4 messages Options
Reply | Threaded
Open this post in threaded view
|

[slicer-devl]Tensor scalar calculation in Python

David Qixiang Chen
Hi,
I've been trying to build a module to output the primary, secondary
and tertiary eigen values as a scalar volume with Python, so I can do
some additional scalar measures such as radial diffusivity.
Following the guide at
http://www.slicer.org/slicerWiki/index.php/Slicer3:Python#Accessing_Tensor_Data,
I was able to get what I believe is the D matrix.
My question is how do I get the 3 eigen values from there.

I looked at the source code at
http://www.na-mic.org/ViewVC/index.cgi/trunk/Libs/vtkTeem/vtkDiffusionTensorMathematics.cxx?view=markup
under TeemEigenSolver(), but it seems to be calling some functions
that I'm not sure where to find in Python.

Thanks!
_______________________________________________
slicer-devel mailing list
[hidden email]
http://massmail.spl.harvard.edu/mailman/listinfo/slicer-devel
To unsubscribe: send email to [hidden email] with unsubscribe as the subject
Reply | Threaded
Open this post in threaded view
|

Re: [slicer-devl]Tensor scalar calculation in Python

pieper
Administrator
Hi David -

You should be able to use the teem.py package bundled with slicer - just
do 'import teem' and then access teem.tenEigensolve_d().  I'm not sure
the right arguments off the top of my head, but it should work.  See
http://teem.sourceforge.net for more info on the teem api.  Let us know
how it goes...

-Steve

David Qixiang Chen wrote:

> Hi,
> I've been trying to build a module to output the primary, secondary
> and tertiary eigen values as a scalar volume with Python, so I can do
> some additional scalar measures such as radial diffusivity.
> Following the guide at
> http://www.slicer.org/slicerWiki/index.php/Slicer3:Python#Accessing_Tensor_Data,
> I was able to get what I believe is the D matrix.
> My question is how do I get the 3 eigen values from there.
>
> I looked at the source code at
> http://www.na-mic.org/ViewVC/index.cgi/trunk/Libs/vtkTeem/vtkDiffusionTensorMathematics.cxx?view=markup
> under TeemEigenSolver(), but it seems to be calling some functions
> that I'm not sure where to find in Python.
>
> Thanks!
> _______________________________________________
> slicer-devel mailing list
> [hidden email]
> http://massmail.spl.harvard.edu/mailman/listinfo/slicer-devel
> To unsubscribe: send email to [hidden email] with unsubscribe as the subject
_______________________________________________
slicer-devel mailing list
[hidden email]
http://massmail.spl.harvard.edu/mailman/listinfo/slicer-devel
To unsubscribe: send email to [hidden email] with unsubscribe as the subject
Reply | Threaded
Open this post in threaded view
|

Re: [slicer-devl]Tensor scalar calculation in Python

Demian Wassermann
In reply to this post by David Qixiang Chen
Hi David,


There is a way of doing it which might be a bit slower than using teem:

Let's suppose that you have the tensors after this line and a voxel  
list, meaning a tuple of three numpy arrays with equal length N and  
type int

voxels = (
        [i0,i1......iN-1],
        [j0,j1.....jN-1],
        [k0,k1....kN-1]
)

tensors = i.reshape(d+[3,3])

if you do

import from numpy import linalg

eigenDecomp = [ linalg.eigh( t ) for v in tensors[ tuple(voxels) ]


eigenDecomp is now a list of tuples where the first element of the  
tuple is the eigenvalue list and the second a matrix where the  
corresponding eigenvectos are the columns.

If you want to have a quicker better representation at the expense of  
writing a bit more:

eigenVals = zeros( (N,3) )
eigenVecs = zeros( (N,3,3) )
selectedTensors = tensors[ voxels ]

for i in xrange(N):
        eigenVals[i],eigenVecs[i] = numpy.eigh( selectedTensors[i] )
       

Hope it helps

D




--
Demian Wassermann
[hidden email]
PhD Student
Odyssee Research Project
INRIA Sophia-Antipolis
2004 route des lucioles - FR-06902




On Dec 14, 2009, at 6:55 PM, David Qixiang Chen wrote:

> Hi,
> I've been trying to build a module to output the primary, secondary
> and tertiary eigen values as a scalar volume with Python, so I can do
> some additional scalar measures such as radial diffusivity.
> Following the guide at
> http://www.slicer.org/slicerWiki/index.php/ 
> Slicer3:Python#Accessing_Tensor_Data,
> I was able to get what I believe is the D matrix.
> My question is how do I get the 3 eigen values from there.
>
> I looked at the source code at
> http://www.na-mic.org/ViewVC/index.cgi/trunk/Libs/vtkTeem/ 
> vtkDiffusionTensorMathematics.cxx?view=markup
> under TeemEigenSolver(), but it seems to be calling some functions
> that I'm not sure where to find in Python.
>
> Thanks!
> _______________________________________________
> slicer-devel mailing list
> [hidden email]
> http://massmail.spl.harvard.edu/mailman/listinfo/slicer-devel
> To unsubscribe: send email to slicer-devel-
> [hidden email] with unsubscribe as the subject

_______________________________________________
slicer-devel mailing list
[hidden email]
http://massmail.spl.harvard.edu/mailman/listinfo/slicer-devel
To unsubscribe: send email to [hidden email] with unsubscribe as the subject
Reply | Threaded
Open this post in threaded view
|

Re: [slicer-devl]Tensor scalar calculation in Python

David Qixiang Chen
Hi,
It took some trials, but I got things to work with teem.tenEigensolve_d()
Thanks for your help!

Here's my code:

        if outputType == 'E1':
                output_array_cmd = 'eigens[0]'
        elif outputType == 'E2':
                output_array_cmd = 'eigens[1]'
        elif outputType == 'E3':
                output_array_cmd = 'eigens[2]'
        elif outputType == 'RD':
                output_array_cmd = '(eigens[1]+eigens[2])/2'

        for i in range(s[0]):
                for j in range(s[1]):
                        for k in range(s[2]):
                                m = tensors[i,j,k]
                                t[1] = m[0][0];
      t[2] = m[0][1];
    t[3] = m[0][2];
  t[4] = m[1][1];
      t[5] = m[1][2];
    t[6] = m[2][2];
                                t_ctypes = ctypes.cast(t.ctypes.data, ctypes.POINTER(ctypes.c_double))

                                teem.tenEigensolve_d(eigens_ctypes, vecs_ctypes, t_ctypes)
                               
                                output_array[i,j,k] = eval(output_array_cmd)


- David


On Mon, Dec 14, 2009 at 5:33 PM, Demian Wassermann
<[hidden email]> wrote:

> Hi David,
>
>
> There is a way of doing it which might be a bit slower than using teem:
>
> Let's suppose that you have the tensors after this line and a voxel list,
> meaning a tuple of three numpy arrays with equal length N and type int
>
> voxels = (
>        [i0,i1......iN-1],
>        [j0,j1.....jN-1],
>        [k0,k1....kN-1]
> )
>
> tensors = i.reshape(d+[3,3])
>
> if you do
>
> import from numpy import linalg
>
> eigenDecomp = [ linalg.eigh( t ) for v in tensors[ tuple(voxels) ]
>
>
> eigenDecomp is now a list of tuples where the first element of the tuple is
> the eigenvalue list and the second a matrix where the corresponding
> eigenvectos are the columns.
>
> If you want to have a quicker better representation at the expense of
> writing a bit more:
>
> eigenVals = zeros( (N,3) )
> eigenVecs = zeros( (N,3,3) )
> selectedTensors = tensors[ voxels ]
>
> for i in xrange(N):
>        eigenVals[i],eigenVecs[i] = numpy.eigh( selectedTensors[i] )
>
>
> Hope it helps
>
> D
>
>
>
>
> --
> Demian Wassermann
> [hidden email]
> PhD Student
> Odyssee Research Project
> INRIA Sophia-Antipolis
> 2004 route des lucioles - FR-06902
>
>
>
>
> On Dec 14, 2009, at 6:55 PM, David Qixiang Chen wrote:
>
>> Hi,
>> I've been trying to build a module to output the primary, secondary
>> and tertiary eigen values as a scalar volume with Python, so I can do
>> some additional scalar measures such as radial diffusivity.
>> Following the guide at
>>
>> http://www.slicer.org/slicerWiki/index.php/Slicer3:Python#Accessing_Tensor_Data,
>> I was able to get what I believe is the D matrix.
>> My question is how do I get the 3 eigen values from there.
>>
>> I looked at the source code at
>>
>> http://www.na-mic.org/ViewVC/index.cgi/trunk/Libs/vtkTeem/vtkDiffusionTensorMathematics.cxx?view=markup
>> under TeemEigenSolver(), but it seems to be calling some functions
>> that I'm not sure where to find in Python.
>>
>> Thanks!
>> _______________________________________________
>> slicer-devel mailing list
>> [hidden email]
>> http://massmail.spl.harvard.edu/mailman/listinfo/slicer-devel
>> To unsubscribe: send email to
>> [hidden email] with unsubscribe as the
>> subject
>
>
_______________________________________________
slicer-devel mailing list
[hidden email]
http://massmail.spl.harvard.edu/mailman/listinfo/slicer-devel
To unsubscribe: send email to [hidden email] with unsubscribe as the subject