dipy logo

Site Navigation

NIPY Community

Table Of Contents

Previous topic

dipy.viz.fvtk

dipy.viz.fvtk

Fvtk module implements simple visualization functions using VTK. Fos means light in Greek.

The main idea is the following: A window can have one or more renderers. A renderer can have none, one or more actors. Examples of actors are a sphere, line, point etc. You basically add actors in a renderer and in that way you can visualize the forementioned objects e.g. sphere, line ...

Examples

>>> from dipy.viz import fvtk
>>> r=fvtk.ren()
>>> a=fvtk.axes()
>>> fvtk.add(r,a)
>>> #fvtk.show(r)
dipy.viz.fvtk.add(ren, a)

Add a specific actor

dipy.viz.fvtk.annotatePick(object, event)

Create a Python function to create the text for the text mapper used to display the results of picking.

dipy.viz.fvtk.axes(scale=(1, 1, 1), colorx=(1, 0, 0), colory=(0, 1, 0), colorz=(0, 0, 1), opacity=1)

Create an actor with the coordinate system axes where red = x, green = y, blue =z.

dipy.viz.fvtk.clear(ren)

Remove all actors from the renderer

dipy.viz.fvtk.colors(v, colormap, auto=True)

Create colors from a specific colormap and return it as an array of shape (N,3) where every row gives the corresponding r,g,b value. The colormaps we use are similar with that of pylab.

Current options for colormaps are ‘jet’,’blues’,’blue_red’, ‘accent’

Notes

If you want to add more colormaps here is what you could do. Go to this website http://www.scipy.org/Cookbook/Matplotlib/Show_colormaps see which colormap you need and then get in pylab using the cm.datad dictionary.

e.g. cm.datad[‘jet’]

{‘blue’: ((0.0, 0.5, 0.5),
(0.11, 1, 1), (0.34000000000000002, 1, 1), (0.65000000000000002, 0, 0), (1, 0, 0)),
‘green’: ((0.0, 0, 0),
(0.125, 0, 0), (0.375, 1, 1), (0.64000000000000001, 1, 1), (0.91000000000000003, 0, 0), (1, 0, 0)),
‘red’: ((0.0, 0, 0),
(0.34999999999999998, 0, 0), (0.66000000000000003, 1, 1), (0.89000000000000001, 1, 1), (1, 0.5, 0.5))}
dipy.viz.fvtk.contour(vol, voxsz=(1.0, 1.0, 1.0), affine=None, levels=[50], colors=[array([ 1., 0., 0.])], opacities=[0.5])

Take a volume and draw surface contours for any any number of thresholds (levels) where every contour has its own color and opacity

Parameters :

vol : array, shape (N, M, K)

an array representing the volumetric dataset for which we will draw some beautiful contours .

voxsz : sequence of 3 floats

default (1., 1., 1.)

affine : not used here

levels : sequence of thresholds for the contours taken from image values

needs to be same datatype as vol

colors : array, shape (N,3) with the rgb values in where r,g,b belong to [0,1]

opacities : sequence of floats [0,1]

Returns :

ass: assembly of actors :

representing the contour surfaces

Examples

>>> import numpy as np
>>> from dipy.viz import fvtk
>>> A=np.zeros((10,10,10))
>>> A[3:-3,3:-3,3:-3]=1
>>> r=fvtk.ren()
>>> fvtk.add(r,fvtk.contour(A,levels=[1]))
>>> #fvtk.show(r)
dipy.viz.fvtk.crossing(a, ind, sph, scale, orient=False)

visualize a volume of crossings

Examples

See ‘dipy/doc/examples/visualize_crossings.py’ at Examples

dipy.viz.fvtk.dots(points, color=(1, 0, 0), opacity=1)

Create one or more 3d dots(points) returns one actor handling all the points

ellipsoid(R=array([[2, 0, 0],
[0, 1, 0],
[0, 0, 1]]), position=(0, 0, 0), thetares=20, phires=20, color=(0, 0, 1), opacity=1, tessel=0)

Create a ellipsoid actor. Stretch a unit sphere to make it an ellipsoid under a 3x3 translation matrix R

R=sp.array([[2, 0, 0],
[0, 1, 0], [0, 0, 1] ])
dipy.viz.fvtk.label(ren, text='Origin', pos=(0, 0, 0), scale=(0.2, 0.2, 0.2), color=(1, 1, 1))

Create a label actor This actor will always face the camera

Parameters :

ren : vtkRenderer() object as returned from ren()

text : a text for the label

pos : left down position of the label

scale : change the size of the label

color : (r,g,b) and RGB tuple

Returns :

vtkActor object :

Examples

>>> from dipy.viz import fvtk  
>>> r=fvtk.ren()    
>>> l=fvtk.label(r)
>>> fvtk.add(r,l)
>>> #fvtk.show(r)
dipy.viz.fvtk.line(lines, colors, opacity=1, linewidth=1)

Create an actor for one or more lines.

Parameters :

lines : list of arrays representing lines as 3d points for example

lines=[np.random.rand(10,3),np.random.rand(20,3)] represents 2 lines the first with 10 points and the second with 20 points in x,y,z coordinates.

colors : array, shape (N,3)

Colormap where every triplet is encoding red, green and blue e.g. r1,g1,b1 r2,g2,b2 ... rN,gN,bN

where 0=<r<=1, 0=<g<=1, 0=<b<=1

opacity : float, default 1

0<=transparency <=1

linewidth : float, default is 1

line thickness

Returns :

vtkActor object :

Examples

>>> from dipy.viz import fvtk
>>> r=fvtk.ren()    
>>> lines=[np.random.rand(10,3),np.random.rand(20,3)]    
>>> colors=np.random.rand(2,3)
>>> c=fvtk.line(lines,colors)    
>>> fvtk.add(r,c)
>>> #fvtk.show(r)
dipy.viz.fvtk.record(ren=None, cam_pos=None, cam_focal=None, cam_view=None, out_path=None, n_frames=10, az_ang=10, magnification=1, size=(300, 300), bgr_color=(0, 0, 0))

This will record a video of your scene

Records a video as a series of .png files of your scene by rotating the azimuth angle az_angle in every frame.

Parameters :

ren : vtkRenderer() object

as returned from function ren()

cam_pos : None or sequence (3,), optional

camera position

cam_focal : None or sequence (3,), optional

camera focal point

cam_view : None or sequence (3,), optional

camera view up

out_path : str, optional

output directory for the frames

n_frames : int, optional

number of frames to save, default 10

az_ang : float, optional

azimuthal angle of camera rotation.

magnification : int, optional

how much to magnify the saved frame

Examples

>>> from dipy.viz import fvtk
>>> r=fvtk.ren()
>>> a=fvtk.axes()    
>>> from dipy.viz import fvtk
>>> r=fvtk.ren()
>>> fvtk.add(r,fvtk.axes())
>>> #uncomment below to record
>>> #fvtk.record(r,cam_pos=(0,0,-10))
dipy.viz.fvtk.ren()

Create a renderer

Returns :a vtkRenderer() object :

Examples

>>> from dipy.viz import fvtk
>>> import numpy as np
>>> r=fvtk.ren()    
>>> lines=[np.random.rand(10,3)]        
>>> c=fvtk.line(lines,fvtk.red)    
>>> fvtk.add(r,c)
>>> #fvtk.show(r)    
dipy.viz.fvtk.rm(ren, a)

Remove a specific actor

dipy.viz.fvtk.rm_all(ren)

Remove all actors from the renderer

dipy.viz.fvtk.show(ren, title='dipy.viz.fvtk', size=(300, 300), png_magnify=1)

Show window

Parameters :

ren : vtkRenderer() object

as returned from function ren()

title : string

a string for the window title bar

size : (int, int)

(width,height) of the window

png_magnify : int

number of times to magnify the screenshot

Notes

If you want to:

  • navigate in the the 3d world use the left - middle - right mouse buttons
  • reset the screen press ‘r’
  • save a screenshot press ‘s’
  • quit press ‘q’

Examples

>>> import numpy as np
>>> from dipy.viz import fvtk
>>> r=fvtk.ren()    
>>> lines=[np.random.rand(10,3),np.random.rand(20,3)]    
>>> colors=np.array([[0.2,0.2,0.2],[0.8,0.8,0.8]])
>>> c=fvtk.line(lines,colors)
>>> fvtk.add(r,c)
>>> l=fvtk.label(r)
>>> fvtk.add(r,l)
>>> #fvtk.show(r)
dipy.viz.fvtk.slicer(ren, vol, voxsz=(1.0, 1.0, 1.0), affine=None, contours=1, planes=1, levels=[20, 30, 40], opacities=[0.8, 0.7, 0.3], colors=None, planesx=[20, 30], planesy=[30, 40], planesz=[20, 30])

Slicer and contour rendering of 3d volumes

Parameters :

vol : array, shape (N, M, K), dtype uint8

an array representing the volumetric dataset that we want to visualize using volumetric rendering

voxsz : sequence of 3 floats

default (1., 1., 1.)

affine : array, shape (4,4), default None

as given by volumeimages

contours : bool 1 to show contours

planes : boolean 1 show planes

levels : contour levels

opacities : opacity for every contour level

colors : None or

planesx : saggital

planesy : coronal

planesz : axial

Examples

>>> import numpy as np
>>> from dipy.viz import fvtk
>>> x, y, z = np.ogrid[-10:10:80j, -10:10:80j, -10:10:80j]
>>> s = np.sin(x*y*z)/(x*y*z)
>>> r=fvtk.ren()    
>>> #fvtk.slicer(r,s) #does showing too 
dipy.viz.fvtk.sphere(position=(0, 0, 0), radius=0.5, thetares=8, phires=8, color=(0, 0, 1), opacity=1, tessel=0)

Create a sphere actor

dipy.viz.fvtk.tube(point1=(0, 0, 0), point2=(1, 0, 0), color=(1, 0, 0), opacity=1, radius=0.1, capson=1, specular=1, sides=8)

Deprecated Wrap a tube around a line connecting point1 with point2 with a specific radius

dipy.viz.fvtk.volume(vol, voxsz=(1.0, 1.0, 1.0), affine=None, center_origin=1, info=0, maptype=0, trilinear=1, iso=0, iso_thr=100, opacitymap=None, colormap=None)

Create a volume and return a volumetric actor using volumetric rendering. This function has many different interesting capabilities. The maptype, opacitymap and colormap are the most crucial parameters here.

Parameters :

vol : array, shape (N, M, K), dtype uint8

an array representing the volumetric dataset that we want to visualize using volumetric rendering

voxsz : sequence of 3 floats

default (1., 1., 1.)

affine : array, shape (4,4), default None

as given by volumeimages

center_origin : int {0,1}, default 1

it considers that the center of the volume is the

point (-vol.shape[0]/2.0+0.5,-vol.shape[1]/2.0+0.5,-vol.shape[2]/2.0+0.5)

info : int {0,1}, default 1

if 1 it prints out some info about the volume, the method and the dataset.

trilinear: int {0,1}, default 1 :

Use trilinear interpolation, default 1, gives smoother rendering. If you want faster interpolation use 0 (Nearest).

maptype : int {0,1}, default 0,

The maptype is a very important parameter which affects the raycasting algorithm in use for the rendering. The options are: If 0 then vtkVolumeTextureMapper2D is used. If 1 then vtkVolumeRayCastFunction is used.

iso : int {0,1} default 0,

If iso is 1 and maptype is 1 then we use vtkVolumeRayCastIsosurfaceFunction which generates an isosurface at the predefined iso_thr value. If iso is 0 and maptype is 1 vtkVolumeRayCastCompositeFunction is used.

iso_thr : int, default 100,

if iso is 1 then then this threshold in the volume defines the value which will be used to create the isosurface.

opacitymap : array, shape (N,2), default None.

The opacity map assigns a transparency coefficient to every point in the volume. The default value uses the histogram of the volume to calculate the opacitymap.

colormap : array, shape (N,4), default None.

The color map assigns a color value to every point in the volume. When None from the histogram it uses a red-blue colormap.

Returns :

vtkVolume :

Notes

What is the difference between TextureMapper2D and RayCastFunction? Coming soon... See VTK user’s guide [book] & The Visualization Toolkit [book] and VTK’s online documentation & online docs.

What is the difference between RayCastIsosurfaceFunction and RayCastCompositeFunction? Coming soon... See VTK user’s guide [book] & The Visualization Toolkit [book] and VTK’s online documentation & online docs.

What about trilinear interpolation? Coming soon... well when time permits really ... :-)

Examples

First example random points

>>> from dipy.viz import fvtk
>>> import numpy as np
>>> vol=100*np.random.rand(100,100,100)
>>> vol=vol.astype('uint8')
>>> print vol.min(), vol.max()
0 99
>>> r = fvtk.ren()
>>> v = fvtk.volume(vol)
>>> fvtk.add(r,v)
>>> #fvtk.show(r)

Second example with a more complicated function

>>> from dipy.viz import fvtk
>>> import numpy as np
>>> x, y, z = np.ogrid[-10:10:20j, -10:10:20j, -10:10:20j]
>>> s = np.sin(x*y*z)/(x*y*z)
>>> r = fvtk.ren()
>>> v = fvtk.volume(s)
>>> fvtk.add(r,v)
>>> #fvtk.show(r)

If you find this function too complicated you can always use mayavi. Please do not forget to use the -wthread switch in ipython if you are running mayavi.

from enthought.mayavi import mlab import numpy as np x, y, z = np.ogrid[-10:10:20j, -10:10:20j, -10:10:20j] s = np.sin(x*y*z)/(x*y*z) mlab.pipeline.volume(mlab.pipeline.scalar_field(s)) mlab.show()

More mayavi demos are available here:

http://code.enthought.com/projects/mayavi/docs/development/html/mayavi/mlab.html