dipy logo

Site Navigation

NIPY Community

Table Of Contents

Previous topic

dipy.viz.fvtk

Next topic

A quick overview of features

dipy.viz.fvtk

Fvtk module implements simple visualization functions using VTK.

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.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 : (N, M, K) ndarray

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

voxsz : (3,) array_like

Voxel size.

affine : None

Not used.

levels : array_like

Sequence of thresholds for the contours taken from image values needs to be same datatype as vol.

colors : (N, 3) ndarray

RGB values in [0,1].

opacities : array_like

Opacities of contours.

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.create_colormap(v, name='jet', 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 those of pylab.

Parameters :

v : (N,) array

vector of values to be mapped in RGB colors according to colormap

name : str. ‘jet’, ‘blues’, ‘blue_red’, ‘accent’

name of the colourmap

auto : bool,

if auto is True then v is interpolated to [0, 10] from v.min() to v.max()

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.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

Renderer as returned by ren().

text : str

Text for the label.

pos : (3,) array_like, optional

Left down position of the label.

scale : (3,) array_like

Changes the size of the label.

color : (3,) array_like

Label color as (r,g,b) tuple.

Returns :

l : vtkActor object

Label.

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, optional

0 <= transparency <= 1

linewidth : float, optional

Line thickness.

Returns :

v : vtkActor object

Line.

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, path_numbering=False, n_frames=10, az_ang=10, magnification=1, size=(300, 300), bgr_color=(0, 0, 0), verbose=False)

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

path_numbering : bool

when recording it changes out_path ot out_path + str(frame number)

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()
>>> fvtk.add(r,a)
>>> #uncomment below to record
>>> #fvtk.record(r)
>>> #check for new images in current directory
dipy.viz.fvtk.ren()

Create a renderer.

Returns :

v : vtkRenderer() object

Renderer.

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

Voxel size.

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

As given by volumeimages.

contours : bool 1 to show contours

Whether to show contours.

planes : boolean 1 show planes

Whether to show planes.

levels : sequence

Contour levels.

opacities : sequence

Opacity for every contour level.

colors : None or sequence of 3-tuples

Color for each contour level.

planesx : (2,) array_like

Saggital.

planesy : (2,) array_like

Coronal.

planesz : :

Axial (2,) array_like

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.sphere_funcs(sphere_values, sphere, image=None, colormap='jet', scale=2.2, norm=True, radial_scale=True)

Plot many morphed spheres simultaneously.

Parameters :

sphere_values : (M,) or (X, M) or (X, Y, M) or (X, Y, Z, M) ndarray

Values on the sphere.

sphere : Sphere

image : None,

Not yet supported.

colormap : None or ‘jet’

If None then no color is used.

scale : float,

Distance between spheres.

norm : bool,

Normalize sphere_values.

radial_scale : bool,

Scale sphere points according to odf values.

Returns :

actor : vtkActor

Spheres.

Examples

>>> from dipy.viz import fvtk
>>> r = fvtk.ren()
>>> odfs = np.ones((5, 5, 724))
>>> odfs[..., 0] = 2.
>>> from dipy.data import get_sphere
>>> sphere = get_sphere('symmetric724')
>>> fvtk.add(r, fvtk.sphere_funcs(odfs, sphere))
>>> #fvtk.show(r)
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 : (3,) array_like

Voxel size.

affine : (4, 4) ndarray

As given by volumeimages.

center_origin : int {0,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}

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

trilinear : int {0,1}

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

maptype : int {0,1}

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}

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

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

opacitymap : (2, 2) ndarray

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 : (4, 4) ndarray

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 :

v : vtkVolume

Volume.

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