librheolef1 6.7-6 / usr / share / rheolef / mayavi2_rheolef.py

##
## This file is part of Rheolef.
##
## Copyright (C) 2000-2009 Pierre Saramito 
##
## Rheolef is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 2 of the License, or
## (at your option) any later version.
##
## Rheolef is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with Rheolef; if not, write to the Free Software
## Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
##
# --------------------------------------------------------------------------
# 2D interaction : vtkInteractorStyleImage
# --------------------------------------------------------------------------
# The bindings keep the camera's view plane normal perpendicular to the x-y plane.
# In summary the mouse events are as follows: 
#	+ Left Mouse button triggers window level events 
#	+ CTRL Left Mouse spins the camera around its view plane normal
#	+ SHIFT Left Mouse pans the camera 
#	+ CTRL SHIFT Left Mouse dollys (a positional zoom) the camera
#	+ Middle mouse button pans the camera
#	+ Right mouse button dollys the camera.
#	+ SHIFT Right Mouse triggers pick events
# Note that the renderer's actors are not moved; instead the camera is moved.
# --------------------------------------------------------------------------
# 3D interaction : vtkInteractorStyleSwitch 
# --------------------------------------------------------------------------
# allows handles interactively switching between four interactor styles:
#	joystick actor
#	joystick camera
#	trackball actor
#	trackball camera.
# Type 'j' or 't' to select joystick or trackball, and type 'c' or 'a' to select camera or actor.
# The default interactor style is joystick camera.
#
# 1) joystick actor
#    manipulate objects in the scene independently of one another
#    allows the user to interact with (rotate, zoom, etc.) separate objects in the scene independent
#    of each other. The position of the mouse relative to the center of the object determines the speed
#    of the object's motion. The mouse's velocity detemines the acceleration of the object's motion, 
#    so the object will continue moving even when the mouse is not moving. For a 3-button mouse,
#    the left button is for rotation, the right button for zooming, the middle button for panning,
#    and ctrl + left button for spinning. (With fewer mouse buttons, ctrl + shift + left button is
#    for zooming, and shift + left button is for panning.)
#
# 2) joystick camera
#    interactive manipulation of the camera
#    allows the user to move (rotate, pan, etc.) the camera, the point of view for the scene.
#    The position of the mouse relative to the center of the scene determines the speed at which
#    the camera moves, and the speed of the mouse movement determines the acceleration of the camera,
#    so the camera continues to move even if the mouse if not moving. For a 3-button mouse, the left
#    button is for rotation, the right button for zooming, the middle button for panning,
#    and ctrl + left button for spinning. (With fewer mouse buttons, ctrl + shift + left button is
#    for zooming, and shift + left button is for panning.)
#
# 3) trackball actor
#    manipulate objects in the scene independent of each other
#    allows the user to interact with (rotate, pan, etc.) objects in the scene indendent of each other.
#    In trackball interaction, the magnitude of the mouse motion is proportional to the actor motion
#    associated with a particular mouse binding. For example, small left-button motions cause small*
#    changes in the rotation of the actor around its center point.
#    The mouse bindings are as follows. For a 3-button mouse, the left button is for rotation,
#    the right button for zooming, the middle button for panning, and ctrl + left button for spinning.
#    (With fewer mouse buttons, ctrl + shift + left button is for zooming,
#    and shift + left button is for panning.)
#
# 4) trackball camera.
#    allows the user to interactively manipulate (rotate, pan, etc.) the camera, the viewpoint of the scene.
#    In trackball interaction, the magnitude of the mouse motion is proportional to the camera motion
#    associated with a particular mouse binding. For example, small left-button motions cause small
#    changes in the rotation of the camera around its focal point.
#    For a 3-button mouse, the left button is for rotation, the right button for zooming,
#    the middle button for panning, and ctrl + left button for spinning.
#    (With fewer mouse buttons, ctrl + shift + left button is for zooming,
#    and shift + left button is for panning.)
#
# --------------------------------------------------------------------------
# includes
# --------------------------------------------------------------------------
from os.path import join, abspath
import numpy	# sqrt, ect
from sys import version_info # mayavi2-4.0 has changed his interface from python-2.7.2...
if version_info[0] >= 2 and version_info[1] >= 7 and version_info[2] >= 2:
    from mayavi.sources.vtk_file_reader import VTKFileReader
    from mayavi.modules.outline import Outline
    from mayavi.modules.axes import Axes
    from mayavi.modules.text import Text
    from mayavi.modules.surface import Surface
    from mayavi.modules.iso_surface import IsoSurface
    from mayavi.modules.scalar_cut_plane import ScalarCutPlane
    from mayavi.modules.volume import Volume
    from mayavi.modules.vectors import Vectors
    from mayavi.modules.tensor_glyph import TensorGlyph
    from mayavi.filters.extract_edges import ExtractEdges
    from mayavi.filters.user_defined import UserDefined
    from mayavi.filters.image_data_probe import ImageDataProbe
    from mayavi.filters.warp_scalar import WarpScalar
    from mayavi.filters.warp_vector import WarpVector
    from mayavi.filters.cut_plane import CutPlane
    from mayavi.components.implicit_plane import ImplicitPlane
    from tvtk.api import tvtk
    from traits.api import Instance, Property, Str, TraitError
else:
    # the same with the old "enthought" prefix (I dont known how to merge these cases)
    from enthought.mayavi.sources.vtk_file_reader import VTKFileReader
    from enthought.mayavi.modules.outline import Outline
    from enthought.mayavi.modules.axes import Axes
    from enthought.mayavi.modules.text import Text
    from enthought.mayavi.modules.surface import Surface
    from enthought.mayavi.modules.iso_surface import IsoSurface
    from enthought.mayavi.modules.scalar_cut_plane import ScalarCutPlane
    from enthought.mayavi.modules.volume import Volume
    from enthought.mayavi.modules.vectors import Vectors
    from enthought.mayavi.modules.tensor_glyph import TensorGlyph
    from enthought.mayavi.filters.extract_edges import ExtractEdges
    from enthought.mayavi.filters.user_defined import UserDefined
    from enthought.mayavi.filters.image_data_probe import ImageDataProbe
    from enthought.mayavi.filters.warp_scalar import WarpScalar
    from enthought.mayavi.filters.warp_vector import WarpVector
    from enthought.mayavi.filters.cut_plane import CutPlane
    from enthought.mayavi.components.implicit_plane import ImplicitPlane
    from enthought.tvtk.api import tvtk
    from enthought.traits.api import Instance, Property, Str, TraitError

# --------------------------------------------------------------------------
# utilities
# --------------------------------------------------------------------------
def start_new_scene (view, option):

    view.engine.new_scene()
    scene = view.engine.current_scene.scene
    scene.disable_render = True

    # bg=white, fg=black
    scene.background = (1,1,1)
    scene.foreground = (0,0,0)

    return scene


def show_scene (scene):

    scene.disable_render = False

def mayavi2_axis (view, option):

    if option['view_2d'] or option['view_1d']:
        option['stereo'] = 0
        scene = view.engine.current_scene.scene
        #print "scene.interactor.interactor_style = ", scene.interactor.interactor_style
        scene.interactor.interactor_style = tvtk.InteractorStyleImage()

    # 2D or 3D and not stereo
    # axes are not rendered as stereo and disturbs the outline
    axes = Axes()
    if option['elevation']:
      # scale in Y or Z direction: use ranges
      axes.axes.use_ranges = True
      axes.axes.ranges = option['ranges']
      scene = view.engine.current_scene.scene
      scene.parallel_projection = True

    if option['stereo']:
        axes.axes.visibility = False
    view.add_module(axes)
    if option['view_2d']:
        # stupid bug: clean y-axis and still show z-axis:
        #   axes.axes.z_axis_visibility = False
        # thus, only clean z label
        axes.axes.z_label = ''
        axes.axes.z_axis_visibility = 0
    if option['view_1d']:
        axes.axes.y_label = ''
        #axes.axes.z_label = ''
        axes.axes.y_axis_visibility = 0
        #axes.axes.z_axis_visibility = 0

    if not option['view_2d']:
	# 3D and stereo : put small axes at bottom
        try:
            # small orientation axes: only works with vtk-4.5...
            from enthought.mayavi.modules.orientation_axes import OrientationAxes
        except ImportError:
            pass
        else:
            small_axes = OrientationAxes()
            #small_axes.marker.set_viewport(0.0, 0.8, 0.2, 1.0)
            view.add_module(small_axes)


def check_stereo_axis (view, option):

    scene = view.engine.current_scene.scene

    if option['axis']:
        outline = Outline()
        view.add_module(outline)
        mayavi2_axis (view, option)

    if not option['view_2d']:

        # select initial view:
        scene.isometric_view()
        # stereo allowed (tape '3' in the window): does not work...
        scene.stereo = True

        # enter directly in stereo, as typing '3' in the window:
        if option['stereo']:
            window = scene.render_window
            window.stereo_render = True
            window.stereo_type   = 'anaglyph' ; # or 'red_blue'
            show_window = 0
            if show_window:
              try:
	        # avoid error message when menu exits after the main window is closed:
                window.edit_traits()
              except BadWindow:
                pass

    if option['view_2d'] or option['view_1d']:

        # select initial view:
        scene.z_plus_view()
    
# --------------------------------------------------------------------------
# geo visualization
# --------------------------------------------------------------------------
# global variables (to all domains):
# ----------------------------------
global_geo_color = [					\
	( 0, 0, 0),				\
	( 0, 1, 0),				\
	( 0, 0, 1),				\
	( 1, 1, 0),				\
	( 1, 0, 1),				\
	( 0, 1, 1),				\
	(0.8900, 0.1500, 0.2100	),		\
	(0.6400, 0.5800, 0.5000 ),		\
	(1.0000, 0.3800, 0.0100	) ]

def mayavi2_geo_domain (view, option, filename, dom_name, idx):
  
    color = global_geo_color[idx % len(global_geo_color)]

    data = VTKFileReader()
    data.initialize(filename)
    view.add_source(data)

    if idx == 0:
        output = data.reader.get_output()
        x_min, x_max, y_min, y_max, z_min, z_max = output.bounds

	# title:
        #print "output = ", output
        #print "output.number_of_points = ", output.number_of_points
        #print "output.number_of_cells = ", output.number_of_cells
        n_vertex  = output.number_of_points
        n_element = output.number_of_cells
	title = "%s : %d elements, %d vertices"%(dom_name, n_element, n_vertex)
        title_text = Text(text=title)
        title_text.x_position = 0.5
        title_text.y_position = 0.9
        view.add_module(title_text)
        title_text.property.font_size = 15
        title_text.property.color = (0,0,0)
        title_text.property.justification = 'centered'
        title_text.property.font_family = 'times'
	vtk_version = float(tvtk.Version().vtk_version[:3])
        if vtk_version > 5.1:
            title_text.actor.set(text_scale_mode='viewport')
        else:
            title_text.actor.set(scaled_text=False)
        title_text.render()

    # creates the cut plane when idx == 0:
    # memorize it in option[] for next calls when idx > 0
    if option['cut']:
        if idx != 0:
	    cut_plane = option['cut-plane']
        else:
	    drag_plane = CutPlane()
	    #print "drag_plane.filters = ", drag_plane.filters
	    ip = drag_plane.filters[0]
            ip.widget.outline_translation = False
            view.add_module(drag_plane)
	    ip._set_normal(option['normal'])
            if option['origin'][0] < 1e+10:
	        ip._set_origin(option['origin'])
            else:
	        ip._set_origin(data.reader.get_output().center)
            ip.render()
            option['cut-plane'] = ip.plane
	    cut_plane = option['cut-plane']

    if idx == 0 and option['full']:
	# NOTE: when cut & full: we do not see edges inside
        edge_filter = ExtractEdges();
        view.add_module(edge_filter)

    if option['cut']:
        polydata_filter = UserDefined(filter=tvtk.GeometryFilter(), name='GeometryFilter')
        view.add_module(polydata_filter)
        clip_filter = UserDefined(filter=tvtk.ClipPolyData(), name='ClipPolyData')
        clip_filter.filter.clip_function = cut_plane
        clip_filter.filter.inside_out = 1
        view.add_module(clip_filter)

    surface = Surface()
    view.add_module(surface)
    #ICI
    surface.actor.mapper.scalar_visibility = False;
    surface.actor.property.color = color
    if option['fill']:
        surface.actor.property.representation = 'surface'
    else:
        surface.actor.property.representation = 'wireframe'

    if option['label']:
        if idx == 0:
          name = "mesh"
        else:
          name = dom_name
        dom_text = Text(text=name)
        dom_text.x_position = 0.01
        dom_text.y_position = 0.01+idx*0.04
        view.add_module(dom_text)
        dom_text.property.font_size = 15
        dom_text.property.color = color
	vtk_version = float(tvtk.Version().vtk_version[:3])
        if vtk_version > 5.1:
            dom_text.actor.set(text_scale_mode='viewport')
        else:
            dom_text.actor.set(scaled_text=False)
        dom_text.render()

    if idx == 0 and option['cut']:
	# re-read file
        data2 = VTKFileReader()
        data2.initialize(filename)
        view.add_source(data2)
        cut_filter = UserDefined(filter=tvtk.Cutter(), name='Cutter')
	cut_filter.filter.cut_function = cut_plane
        view.add_filter(cut_filter)
        show_cut_surf = option['fill']
        show_cut_grid = not option['fill'] or option['full']
	if show_cut_surf:
            cut_surf = Surface()
            cut_surf.actor.property.representation = 'surface'
            #cut_surf.actor.property.color = (0.9,0.9,0.9)
            cut_surf.actor.property.color = (1, 0, 0)
            view.add_module(cut_surf)
	if show_cut_grid:
            cut_grid = Surface()
            cut_grid.actor.property.representation = 'wireframe'
            cut_grid.actor.property.color = (1,0,0)
            view.add_module(cut_grid)

def mayavi2_geo (view, option, d):

    scene = start_new_scene (view, option)

    has_only_bdry = option['fill'] and not (option['full'] or option['cut'])
    if not has_only_bdry:
        mayavi2_geo_domain (view, option, d[0]+".vtk", d[0], 0)
    for i in range(1,len(d)):
        mayavi2_geo_domain (view, option, d[0]+"."+d[i]+".vtk", d[i], i)

    check_stereo_axis (view,option)
    show_scene (scene)

# --------------------------------------------------------------------------
# non-interactive cut: from file to file
# --------------------------------------------------------------------------
def plane_cutter_on_file (input_file_name, output_file_name, origin, normal, tolerance):
    # --------------------
    # 1. reader
    # --------------------
    data = VTKFileReader()
    data.initialize(input_file_name)
    # --------------------
    # 2. plane
    # --------------------
    plane2 = tvtk.Plane()
    if origin[0] < 1e+10:
      plane2.origin = origin
    else:
      plane2.origin = data.reader.get_output().center
    plane2.normal = normal
    #print "plane2=",plane2
    # --------------------
    # 3. cut filter
    # --------------------
    cutter = UserDefined(filter=tvtk.Cutter(), name='Cutter')
    cutter.filter.cut_function = plane2
    cutter.filter.input  = data.reader.get_output()
    # --------------------
    # 4. merge points
    # --------------------
    cleaner = UserDefined(filter=tvtk.CleanPolyData(), name='Cleaner')
    cleaner.filter.input     = cutter.filter.output
    cleaner.filter.tolerance = tolerance
    # --------------------
    # 5. force triangle
    # --------------------
    caster = UserDefined(filter=tvtk.CastToConcrete(), name='Caster')
    caster.filter.input     = cleaner.filter.output
    triangler = UserDefined(filter=tvtk.TriangleFilter(), name='Triangler')
    triangler.filter.input  = caster.filter.output
    # --------------------
    # F. writer
    # --------------------
    writer = tvtk.DataSetWriter()
    writer.input = triangler.filter.output
    writer.file_name = output_file_name
    writer.write()

# call example:
#input_file_name = "contraction-zr-P2.vtk"
#output_file_name = "cut.vtk"
#origin = (0.0, 1e-5, 0.0)
#normal = (0,1,0)
#tolerance = 1e-14
#plane_cutter_on_file(input_file_name, output_file_name, origin, normal, tolerance)

# --------------------------------------------------------------------------
# scalar field visualization
# --------------------------------------------------------------------------

# build an isovalue table
def build_isovalue_table (fmin, fmax, option):

    nn = option['n_isovalue_negative']
    np = option['n_isovalue'] - nn
    values = [fmin]
    for i in range(1, nn):
      v = fmin - (1.0*i)/(1.0*nn)*fmin
      values = values + [v]
    values = values + [0]
    for i in range(1, np):
      v = (1.0*i)/(1.0*np)*fmax
      values = values + [v]
    values = values + [fmax]
    return values

def mayavi2_field_scalar (view, basename, option):

    #print "ici scalaire..."
    scene = start_new_scene (view, option)

    data = VTKFileReader()
    data.initialize(basename + ".vtk")
    view.add_source(data)
    output = data.reader.get_output()
    #print "output = ", output
    # P1 = point_data data:
    #print "output.point_data = ", output.point_data
    #print "output.point_data.number_of_arrays = ", output.point_data.number_of_arrays
    #print "output.point_data.scalars = ", output.point_data.scalars
    #print "output.point_data.scalars.name = ", output.point_data.scalars.name
    # P0 = cell data:
    #print "output.cell_data = ", output.cell_data
    #print "output.cell_data.scalars = ", output.cell_data.scalars
    #print "output.cell_data.number_of_arrays = ", output.cell_data.number_of_arrays
    #print "output.cell_data.scalars.name = ", output.cell_data.scalars.name
    # determine whether continuous or discontinuous data:
    if output.point_data.number_of_arrays > 0:
        option['approx'] = "P1"
        option['name'] = output.point_data.scalars.name
    else:
        option['approx'] = "P0"
        option['name'] = output.cell_data.scalars.name

    x_min, x_max, y_min, y_max, z_min, z_max = output.bounds
    if (z_max - z_min == 0) and (y_max - y_min != 0) and option['elevation'] and option['approx'] == "P0":
        print "mayavi2_rheolef.py: elevation for P0 approximation not yet supported"
        option['elevation'] = 0
    if (y_max - y_min == 0):
        option['elevation'] = 1;
    if ((z_max - z_min == 0) and option['cut']):
      option['elevation'] = 1
    option['view_2d'] = ((z_max - z_min == 0) and not option['elevation'])
    option['view_1d'] = (z_max - z_min == 0 and y_max - y_min == 0)

    if (z_max - z_min == 0 and y_max - y_min == 0):
	# 1D data:
        #print "1D data..."

        elevation = WarpScalar()
        view.add_module(elevation)
        
        surface = Surface()
        view.add_module(surface)
        surface.enable_contours = True;
        surface.contour.filled_contours = True;
	surface.actor.mapper.scalar_visibility = False;

	mm = surface.module_manager
	ms = mm.scalar_lut_manager
        #print "ms=",ms
        ms.show_legend = False
        ms.render()
	range = tuple(ms.data_range)
        f_min = range[0]
        f_max = range[1]
	if (f_max-f_min > 0):
            scale = 0.5*(x_max-x_min)/(f_max - f_min)
	else:
            scale = 1.0
 	scale = option['scale']*scale
        elevation.filter.scale_factor = scale
        elevation.filter.normal = 0, 1, 0; # 2d-elevation: y=f(x)
        elevation.filter.use_normal = True
 	option['ranges'] = (x_min, x_max, f_min, f_max, 0, 0)
        elevation.render()
        scene.z_plus_view()

    elif (z_max - z_min == 0):
	# 2D data
        if option['elevation']:
            elevation = WarpScalar()
            view.add_module(elevation)
        surface = Surface()
        view.add_module(surface)
	mm = surface.module_manager
	range = tuple(mm.scalar_lut_manager.data_range)
        f_min = range[0]
        f_max = range[1]
        if option['approx'] != "P0" and not option['view_map']:
            surface.enable_contours = True;
            surface.contour.number_of_contours = option['n_isovalue'];	
	    if option['n_isovalue_negative'] == 0:	
                surface.contour.auto_contours = True
	    else:
                surface.contour.auto_contours = False
                surface.contour.contours = build_isovalue_table(f_min,f_max,option)
                #print "surface.contour.contours = ", surface.contour.contours
            if (     (option['fill'] or option['elevation'])
             and (option['color'] != 'black_and_white' )):
                surface.contour.filled_contours = True;
        if option['elevation'] and not option['fill']:
            surface.actor.property.representation = 'wireframe'
        if option['color'] == 'black_and_white':
	    surface.actor.mapper.scalar_visibility = False;
	#print "surface = ", surface
	#print "surface.actor.mapper = ", surface.actor.mapper
        if option['color'] == 'gray':
            mm.scalar_lut_manager.lut_mode = 'gray'
        if option['approx'] != "P0" and not option['view_map']:
            mm.scalar_lut_manager.number_of_colors = option['n_isovalue']
            mm.scalar_lut_manager.number_of_labels = option['n_isovalue']+1
        #print "ici 2D..."
        if option['elevation']:
            #print "ici 2D elevation..."
	    if (f_max-f_min > 0):
                scale = 0.5*numpy.sqrt((x_max-x_min)**2+(y_max-y_min)**2)/(f_max - f_min)
	    else:
                scale = 1.0
 	    scale = option['scale']*scale
            elevation.filter.scale_factor = scale
 	    option['ranges'] = (x_min, x_max, y_min, y_max, f_min, f_max)
            elevation.render()
            # add black lines, usefull when stereo rending
            iso_surface = IsoSurface()
	    iso_surface.contour.auto_contours = True
	    iso_surface.contour.number_of_contours = option['n_isovalue']
 	    iso_surface.actor.mapper.scalar_visibility = False
 	    iso_surface.actor.property.line_width = 3
            view.add_module(iso_surface)
        if option['cut']:
            scalar_cut_plane = ScalarCutPlane()
            scalar_cut_plane.implicit_plane.widget.outline_translation = False
            view.add_module(scalar_cut_plane)
	    scalar_cut_plane.implicit_plane._set_normal(option['normal'])
            if option['origin'][0] < 1e+10:
	      scalar_cut_plane.implicit_plane._set_origin(option['origin'])
 	    scalar_cut_plane.actor.mapper.scalar_visibility = False
 	    prop = scalar_cut_plane.actor.property
	    #print "prop=",prop
 	    prop.line_width = 4

    else:
	# 3D
        view.engine.current_scene.scene.isometric_view()
        if option['volume']:
            structured_filter = ImageDataProbe()
            view.add_module(structured_filter)
            volume = Volume()
	    #print "volume.property = ", volume.property
            #print "volume.volume_property = ", volume.volume_property
            view.add_module(volume)
	    mm = volume.module_manager
        else: 
	    mm = None
            if option['view_map']:
                surface = Surface()
	        # avoid smoothed iso-surface contours on 3D surface:
 		surface.actor.mapper.interpolate_scalars_before_mapping = True
                view.add_module(surface)
	        mm = surface.module_manager
                option['cut'] = 0
            if (option['iso'] or option['view_map']) and option['approx'] == "P1":
                iso_surface = IsoSurface()
                if option['view_map']:
		    iso_surface.contour.auto_contours = True
		    iso_surface.contour.number_of_contours = option['n_isovalue']+2
 		    iso_surface.actor.mapper.scalar_visibility = False
 		    iso_surface.actor.property.line_width = 3

                view.add_module(iso_surface)
                if not option['fill']:
                    iso_surface.actor.property.representation = 'wireframe'
	        if mm == None:
                    mm = iso_surface.module_manager
	    if mm != None:
                if option['color'] == 'gray':
                    mm.scalar_lut_manager.lut_mode = 'gray'
                if option['approx'] != "P0":
                    mm.scalar_lut_manager.number_of_colors = option['n_isovalue']
                    mm.scalar_lut_manager.number_of_labels = option['n_isovalue']+1
            if option['cut']:
                scalar_cut_plane = ScalarCutPlane()
                scalar_cut_plane.implicit_plane.widget.outline_translation = False
                view.add_module(scalar_cut_plane)
		scalar_cut_plane.implicit_plane._set_normal(option['normal'])
                if option['origin'][0] < 1e+10:
		    scalar_cut_plane.implicit_plane._set_origin(option['origin'])
                scalar_cut_plane.implicit_plane.render()
                if option['approx'] != "P0":
		    scalar_cut_plane.enable_contours = True
		    scalar_cut_plane.contour.filled_contours = True
		    scalar_cut_plane.contour.auto_contours = True
		    scalar_cut_plane.contour.number_of_contours = option['n_isovalue']
	        if mm == None:
                    mm = scalar_cut_plane.module_manager
            if option['approx'] == "P0":
                if option['cut']:
                    polydata_filter = UserDefined(filter=tvtk.GeometryFilter(), name='GeometryFilter')
                    view.add_module(polydata_filter)
                    clip_filter = UserDefined(filter=tvtk.ClipPolyData(), name='ClipPolyData')
                    cut_plane = scalar_cut_plane.implicit_plane.plane
                    clip_filter.filter.clip_function = cut_plane
                    clip_filter.filter.inside_out = 1
                    view.add_module(clip_filter)
                surface = Surface()
                view.add_module(surface)
	        if mm == None:
	            mm = surface.module_manager

    #print "mm.scalar_lut_manager = ", mm.scalar_lut_manager
    if not option['view_1d'] and option['color'] != 'black_and_white' and mm != None:
        mm.scalar_lut_manager.show_scalar_bar = True
        if option['label'] != "output" and (option['name'] != "scalars" and option['name'] != "scalar" and option['name'] != "unamed"):
            title = option['name'] + ": " + option['label']
        elif option['label'] != "output":
            title = option['label']
        elif (option['name'] != "scalars" and option['name'] != "scalar" and option['name'] != "unamed"):
            title = option['name']
        else:
            title = ""
        mm.scalar_lut_manager.data_name = title

    check_stereo_axis (view,option)
    show_scene (scene)

# --------------------------------------------------------------------------
# vector field visualization
# --------------------------------------------------------------------------
def h_moy (data):
    #print "data.reader.get_output = ", data.reader.get_output()
    output = data.reader.get_output()
    x_min, x_max, y_min, y_max, z_min, z_max = output.bounds
    #n = output.number_of_points
    n = output.number_of_cells
    if z_max-z_min > 0:
      h0 = numpy.sqrt(((x_max-x_min)**2 + (y_max-y_min)**2 +
                           (z_max-z_min)**2)/(4.0*n**(1.0/3.0)))
    else:
        if y_max-y_min > 0:
          h0 = numpy.sqrt(0.5*((x_max-x_min)**2 + (y_max-y_min)**2))/numpy.sqrt(0.5*n)
        else:
          h0 = (x_max-x_min)/(1.0*n)
    return h0

def mayavi2_field_velocity(view, data, basename, option):

    velocity = Vectors()
    view.add_module(velocity)

    h0 = h_moy(data)
    alpha = h0*option['scale']
    velocity.glyph.glyph.scale_factor = alpha

    mm = velocity.module_manager
    mm.scalar_lut_manager.show_scalar_bar = True
    mm.scalar_lut_manager.data_name       = option['label']

def mayavi2_field_deformation(view, data, basename, option):

    deformation_filter = WarpVector()
    view.add_module(deformation_filter)
    deformation_filter.filter.scale_factor = option['scale']

    surface = Surface()
    if not option['fill']:
        surface.actor.property.representation = 'wireframe'
    else:
        surface.actor.property.representation = 'surface'
    view.add_module(surface)

    mm = surface.module_manager
    mm.scalar_lut_manager.show_scalar_bar = True
    mm.scalar_lut_manager.data_name       = option['label']

def mayavi2_field_vector(view, basename, option):

    scene = start_new_scene (view, option)

    data = VTKFileReader()
    data.initialize(basename + ".vtk")
    view.add_source(data)
    output = data.reader.get_output()
    x_min, x_max, y_min, y_max, z_min, z_max = output.bounds
    option['view_2d'] = (z_max - z_min == 0)
    check_stereo_axis (view,option)

    if option['style'] == 'deformation':
        mayavi2_field_deformation(view, data, basename, option)
    else:
        mayavi2_field_velocity(view, data, basename, option)

    show_scene (scene)

# --------------------------------------------------------------------------
# tensor field visualization
# --------------------------------------------------------------------------
def mayavi2_field_tensor(view, basename, option):

    scene = start_new_scene (view, option)

    data = VTKFileReader()
    data.initialize(basename + ".vtk")
    view.add_source(data)
    output = data.reader.get_output()
    x_min, x_max, y_min, y_max, z_min, z_max = output.bounds
    option['view_2d'] = (z_max - z_min == 0)
    check_stereo_axis (view,option)

    tensor = TensorGlyph()
    h0 = h_moy(data)
    alpha = h0*option['scale']
    tensor.glyph.glyph.scale_factor = alpha
    view.add_module(tensor)

    mm = tensor.module_manager
    mm.scalar_lut_manager.show_scalar_bar = True
    mm.scalar_lut_manager.data_name       = option['label']

    show_scene (scene)