Back to index

salome-paravis  6.5.0
presentations.py
Go to the documentation of this file.
00001 # Copyright (C) 2010-2012  CEA/DEN, EDF R&D
00002 #
00003 # This library is free software; you can redistribute it and/or
00004 # modify it under the terms of the GNU Lesser General Public
00005 # License as published by the Free Software Foundation; either
00006 # version 2.1 of the License.
00007 #
00008 # This library is distributed in the hope that it will be useful,
00009 # but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 # Lesser General Public License for more details.
00012 #
00013 # You should have received a copy of the GNU Lesser General Public
00014 # License along with this library; if not, write to the Free Software
00015 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 #
00017 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 #
00019 
00020 """
00021 This module is intended to provide Python API for building presentations
00022 typical for Post-Pro module (Scalar Map, Deformed Shape, Vectors, etc.)
00023 """
00024 
00025 
00026 from __future__ import division
00027 ##from __future__ import print_function
00028 
00029 import os
00030 import re
00031 import warnings
00032 from math import sqrt, sin, cos, radians
00033 from string import upper
00034 
00035 try:
00036     import pvsimple as pv
00037     # TODO(MZN): to be removed (issue with Point Sprite texture)
00038     #import paravisSM as sm
00039 except:
00040     import paraview.simple as pv
00041     import paraview.servermanager as sm
00042 
00043 
00044 # Constants
00045 EPS = 1E-3
00046 FLT_MIN = 1E-37
00047 VTK_LARGE_FLOAT = 1E+38
00048 GAP_COEFFICIENT = 0.0001
00049 
00050 
00051 # Globals
00052 _current_bar = None
00053 
00054 
00055 # Enumerations
00056 class PrsTypeEnum:
00057     """
00058     Post-Pro presentation types.
00059     """
00060     MESH = 0
00061     SCALARMAP = 1
00062     ISOSURFACES = 2
00063     CUTPLANES = 3
00064     CUTLINES = 4
00065     DEFORMEDSHAPE = 5
00066     DEFORMEDSHAPESCALARMAP = 6
00067     VECTORS = 7
00068     PLOT3D = 8
00069     STREAMLINES = 9
00070     GAUSSPOINTS = 10
00071 
00072     _type2name = {MESH: 'Mesh',
00073                   SCALARMAP: 'Scalar Map',
00074                   ISOSURFACES: 'Iso Surfaces',
00075                   CUTPLANES: 'Cut Planes',
00076                   CUTLINES: 'Cut Lines',
00077                   DEFORMEDSHAPE: 'Deformed Shape',
00078                   DEFORMEDSHAPESCALARMAP: 'Deformed Shape And Scalar Map',
00079                   VECTORS: 'Vectors',
00080                   PLOT3D: 'Plot3D',
00081                   STREAMLINES: 'Stream Lines',
00082                   GAUSSPOINTS: 'Gauss Points'}
00083 
00084     @classmethod
00085     def get_name(cls, type):
00086         """Return presentaion name by its type."""
00087         return cls._type2name[type]
00088 
00089 
00090 class EntityType:
00091     """
00092     Entity types.
00093     """
00094     NODE = 0
00095     CELL = 1
00096 
00097     _type2name = {NODE: 'OnPoint',
00098                   CELL: 'OnCell'}
00099 
00100     _name2type = {'OnPoint': NODE,
00101                   'OnCell': CELL}
00102 
00103     _type2pvtype = {NODE: 'POINT_DATA',
00104                     CELL: 'CELL_DATA'}
00105 
00106     @classmethod
00107     def get_name(cls, type):
00108         """Return entity name (used in full group names) by its type."""
00109         return cls._type2name[type]
00110 
00111     @classmethod
00112     def get_type(cls, name):
00113         """Return entity type by its name (used in full group names)."""
00114         return cls._name2type[name]
00115 
00116     @classmethod
00117     def get_pvtype(cls, type):
00118         """Return entity type from ['CELL_DATA', 'POINT_DATA']"""
00119         return cls._type2pvtype[type]
00120 
00121 
00122 class Orientation:
00123     """Orientation types.
00124 
00125     Defines a set of plane orientation possibilities:
00126       AUTO: plane orientation should be calculated.
00127       XY: plane formed by X and Y axis.
00128       YZ: plane formed by Y and Z axis.
00129       ZX: plane formed by Z and X axis
00130 
00131     """
00132     AUTO = 0
00133     XY = 1
00134     YZ = 2
00135     ZX = 3
00136 
00137 
00138 class GlyphPos:
00139     """Glyph positions.
00140 
00141     Set of elements defining the position of the vector head:
00142       CENTER: in the center of the vector
00143       TAIL: in the tail of the vector
00144       HEAD: in the head of the vector
00145 
00146     """
00147     CENTER = 0
00148     TAIL = 1
00149     HEAD = 2
00150 
00151 
00152 class GaussType:
00153     """
00154     Gauss Points primitive types.
00155     """
00156     SPRITE = 0
00157     POINT = 1
00158     SPHERE = 2
00159 
00160     _type2mode = {SPRITE: 'Texture',
00161                   POINT: 'SimplePoint',
00162                   SPHERE: 'Sphere (Texture)'}
00163 
00164     @classmethod
00165     def get_mode(cls, type):
00166         """Return paraview point sprite mode by the primitive type."""
00167         return cls._type2mode[type]
00168 
00169 
00170 # Auxiliary functions
00171 def process_prs_for_test(prs, view, picture_name, show_bar=True):
00172     """Show presentation and record snapshot image.
00173 
00174     Arguments:
00175       prs: the presentation to show
00176       view: the render view
00177       picture_name: the full name of the graphics file to save
00178       show_bar: to show scalar bar or not
00179 
00180     """
00181     # Show the presentation only
00182     display_only(prs, view)
00183 
00184     # Show scalar bar
00185     if show_bar and _current_bar:
00186         _current_bar.Visibility = 1
00187 
00188     # Reset the view
00189     reset_view(view)
00190 
00191     # Create a directory for screenshot if necessary
00192     file_name = re.sub("\s+", "_", picture_name)
00193     pic_dir = os.path.dirname(picture_name)
00194     if not os.path.exists(pic_dir):
00195         os.makedirs(pic_dir)
00196 
00197     # Save picture
00198     pv.WriteImage(file_name, view=view, Magnification=1)
00199 
00200 
00201 def reset_view(view=None):
00202     """Reset the view.
00203 
00204     Set predefined (taken from Post-Pro) camera settings.
00205     If the view is not passed, the active view is used.
00206 
00207     """
00208     if not view:
00209         view = pv.GetRenderView()
00210 
00211     # Camera preferences
00212     view.CameraFocalPoint = [0.0, 0.0, 0.0]
00213     view.CameraViewUp = [0.0, 0.0, 1.0]
00214     view.CameraPosition = [738.946, -738.946, 738.946]
00215 
00216     # Turn on the headligth
00217     view.LightSwitch = 1
00218     view.LightIntensity = 0.5
00219 
00220     # Use parallel projection
00221     view.CameraParallelProjection = 1
00222 
00223     view.ResetCamera()
00224     pv.Render(view=view)
00225 
00226 
00227 def hide_all(view, to_remove=False):
00228     """Hide all representations in the view."""
00229     if not view:
00230         view = pv.GetRenderView()
00231 
00232     rep_list = view.Representations
00233     for rep in rep_list:
00234         if hasattr(rep, 'Visibility') and rep.Visibility != 0:
00235             rep.Visibility = 0
00236         if to_remove:
00237             view.Representations.remove(rep)
00238     pv.Render(view=view)
00239 
00240 
00241 def display_only(prs, view=None):
00242     """Display only the given presentation in the view."""
00243     hide_all(view)
00244     if (hasattr(prs, 'Visibility') and prs.Visibility != 1):
00245         prs.Visibility = 1
00246     pv.Render(view=view)
00247 
00248 
00249 def set_visible_lines(xy_prs, lines):
00250     """Set visible only the given lines for XYChartRepresentation."""
00251     sv = xy_prs.GetProperty("SeriesVisibilityInfo").GetData()
00252     visible = '0'
00253 
00254     for i in xrange(0, len(sv)):
00255         if i % 2 == 0:
00256             line_name = sv[i]
00257             if line_name in lines:
00258                 visible = '1'
00259             else:
00260                 visible = '0'
00261         else:
00262             sv[i] = visible
00263 
00264     xy_prs.SeriesVisibility = sv
00265 
00266 
00267 def check_vector_mode(vector_mode, nb_components):
00268     """Check vector mode.
00269 
00270     Check if vector mode is correct for the data array with the
00271     given number of components.
00272 
00273     Arguments:
00274       vector_mode: 'Magnitude', 'X', 'Y' or 'Z'
00275       nb_components: number of component in the data array
00276 
00277     Raises:
00278       ValueError: in case of the vector mode is unexistent
00279       or nonapplicable.
00280 
00281     """
00282     if vector_mode not in ('Magnitude', 'X', 'Y', 'Z'):
00283         raise ValueError("Unexistent vector mode: " + vector_mode)
00284 
00285     if ((nb_components == 1 and (vector_mode == 'Y' or vector_mode == 'Z')) or
00286         (nb_components == 2 and  vector_mode == 'Z')):
00287         raise ValueError("Incorrect vector mode " + vector_mode + " for " +
00288                          nb_components + "-component field")
00289 
00290 
00291 def get_vector_component(vector_mode):
00292     """Get vector component as ineger.
00293 
00294     Translate vector component notation from string
00295     to integer:
00296       'Magnitude': -1
00297       'X': 0
00298       'Y': 1
00299       'Z': 2
00300 
00301     """
00302     vcomponent = -1
00303 
00304     if vector_mode == 'X':
00305         vcomponent = 0
00306     elif vector_mode == 'Y':
00307         vcomponent = 1
00308     elif vector_mode == 'Z':
00309         vcomponent = 2
00310 
00311     return vcomponent
00312 
00313 
00314 def get_data_range(proxy, entity, field_name, vector_mode='Magnitude',
00315                    cut_off=False):
00316     """Get data range for the field.
00317 
00318     Arguments:
00319       proxy: the pipeline object, containig data array for the field
00320       entity: the field entity
00321       field_name: the field name
00322       vector_mode: the vector mode ('Magnitude', 'X', 'Y' or 'Z')
00323 
00324     Returns:
00325       Data range as [min, max]
00326 
00327     """
00328     entity_data_info = None
00329     field_data = proxy.GetFieldDataInformation()
00330 
00331     if field_name in field_data.keys():
00332         entity_data_info = field_data
00333     elif entity == EntityType.CELL:
00334         entity_data_info = proxy.GetCellDataInformation()
00335     elif entity == EntityType.NODE:
00336         entity_data_info = proxy.GetPointDataInformation()
00337 
00338     data_range = []
00339 
00340     if field_name in entity_data_info.keys():
00341         vcomp = get_vector_component(vector_mode)
00342         data_range = entity_data_info[field_name].GetComponentRange(vcomp)
00343     else:
00344         pv_entity = EntityType.get_pvtype(entity)
00345         warnings.warn("Field " + field_name +
00346                       " is unknown for " + pv_entity + "!")
00347 
00348     # Cut off the range
00349     if cut_off and (data_range[0] <= data_range[1]):
00350         data_range = list(data_range)
00351         delta = abs(data_range[1] - data_range[0]) * GAP_COEFFICIENT
00352         data_range[0] += delta
00353         data_range[1] -= delta
00354 
00355     return data_range
00356 
00357 
00358 def get_bounds(proxy):
00359     """Get bounds of the proxy in 3D."""
00360     dataInfo = proxy.GetDataInformation()
00361     bounds_info = dataInfo.GetBounds()
00362     return bounds_info
00363 
00364 
00365 def get_x_range(proxy):
00366     """Get X range of the proxy bounds in 3D."""
00367     bounds_info = get_bounds(proxy)
00368     return bounds_info[0:2]
00369 
00370 
00371 def get_y_range(proxy):
00372     """Get Y range of the proxy bounds in 3D."""
00373     bounds_info = get_bounds(proxy)
00374     return bounds_info[2:4]
00375 
00376 
00377 def get_z_range(proxy):
00378     """Get Z range of the proxy bounds in 3D."""
00379     bounds_info = get_bounds(proxy)
00380     return bounds_info[4:6]
00381 
00382 
00383 def is_planar_input(proxy):
00384     """Check if the given input is planar."""
00385     bounds_info = get_bounds(proxy)
00386 
00387     if (abs(bounds_info[0] - bounds_info[1]) <= FLT_MIN or
00388         abs(bounds_info[2] - bounds_info[3]) <= FLT_MIN or
00389         abs(bounds_info[4] - bounds_info[5]) <= FLT_MIN):
00390         return True
00391 
00392     return False
00393 
00394 
00395 def is_data_on_cells(proxy, field_name):
00396     """Check the existence of a field on cells with the given name."""
00397     cell_data_info = proxy.GetCellDataInformation()
00398     return (field_name in cell_data_info.keys())
00399 
00400 
00401 def is_empty(proxy):
00402     """Check if the object contains any points or cells.
00403 
00404     Returns:
00405       True: if the given proxy doesn't contain any points or cells
00406       False: otherwise
00407 
00408     """
00409     data_info = proxy.GetDataInformation()
00410 
00411     nb_cells = data_info.GetNumberOfCells()
00412     nb_points = data_info.GetNumberOfPoints()
00413 
00414     return not(nb_cells + nb_points)
00415 
00416 
00417 def get_orientation(proxy):
00418     """Get the optimum cutting plane orientation for Plot 3D."""
00419     orientation = Orientation.XY
00420 
00421     bounds = get_bounds(proxy)
00422     delta = [bounds[1] - bounds[0],
00423              bounds[3] - bounds[2],
00424              bounds[5] - bounds[4]]
00425 
00426     if (delta[0] >= delta[1] and delta[0] >= delta[2]):
00427         if (delta[1] >= delta[2]):
00428             orientation = Orientation.XY
00429         else:
00430             orientation = Orientation.ZX
00431     elif (delta[1] >= delta[0] and delta[1] >= delta[2]):
00432         if (delta[0] >= delta[2]):
00433             orientation = Orientation.XY
00434         else:
00435             orientation = Orientation.YZ
00436     elif (delta[2] >= delta[0] and delta[2] >= delta[1]):
00437         if (delta[0] >= delta[1]):
00438             orientation = Orientation.ZX
00439         else:
00440             orientation = Orientation.YZ
00441 
00442     return orientation
00443 
00444 
00445 def dot_product(a, b):
00446     """Dot product of two 3-vectors."""
00447     dot = a[0] * b[0] + a[1] * b[1] + a[2] * b[2]
00448     return dot
00449 
00450 
00451 def multiply3x3(a, b):
00452     """Mutltiply one 3x3 matrix by another."""
00453     c = [[0, 0, 0],
00454          [0, 0, 0],
00455          [0, 0, 0]]
00456 
00457     for i in xrange(3):
00458         c[0][i] = a[0][0] * b[0][i] + a[0][1] * b[1][i] + a[0][2] * b[2][i]
00459         c[1][i] = a[1][0] * b[0][i] + a[1][1] * b[1][i] + a[1][2] * b[2][i]
00460         c[2][i] = a[2][0] * b[0][i] + a[2][1] * b[1][i] + a[2][2] * b[2][i]
00461 
00462     return c
00463 
00464 
00465 def get_rx(ang):
00466     """Get X rotation matrix by angle."""
00467     rx = [[1.0, 0.0,      0.0],
00468           [0.0, cos(ang), -sin(ang)],
00469           [0.0, sin(ang), cos(ang)]]
00470 
00471     return rx
00472 
00473 
00474 def get_ry(ang):
00475     """Get Y rotation matrix by angle."""
00476     ry = [[cos(ang),  0.0, sin(ang)],
00477           [0.0,       1.0, 0.0],
00478           [-sin(ang), 0.0, cos(ang)]]
00479 
00480     return ry
00481 
00482 
00483 def get_rz(ang):
00484     """Get Z rotation matrix by angle."""
00485     rz = [[cos(ang), -sin(ang), 0.0],
00486           [sin(ang), cos(ang),  0.0],
00487           [0.0,      0.0,       1.0]]
00488 
00489     return rz
00490 
00491 
00492 def get_normal_by_orientation(orientation, ang1=0, ang2=0):
00493     """Get normal for the plane by its orientation."""
00494     i_plane = 0
00495     rotation = [[], [], []]
00496     rx = ry = rz = [[1.0, 0.0, 0.0],
00497                     [0.0, 1.0, 0.0],
00498                     [0.0, 0.0, 1.0]]
00499 
00500     normal = [0.0, 0.0, 0.0]
00501     if orientation == Orientation.XY:
00502         if abs(ang1) > EPS:
00503             rx = get_rx(ang1)
00504         if abs(ang2) > EPS:
00505             ry = get_ry(ang2)
00506         rotation = multiply3x3(rx, ry)
00507         i_plane = 2
00508     elif orientation == Orientation.ZX:
00509         if abs(ang1) > EPS:
00510             rz = get_rz(ang1)
00511         if abs(ang2) > EPS:
00512             rx = get_rx(ang2)
00513         rotation = multiply3x3(rz, rx)
00514         i_plane = 1
00515     elif orientation == Orientation.YZ:
00516         if abs(ang1) > EPS:
00517             ry = get_ry(ang1)
00518         if abs(ang2) > EPS:
00519             rz = get_rz(ang2)
00520         rotation = multiply3x3(ry, rz)
00521         i_plane = 0
00522 
00523     for i in xrange(0, 3):
00524         normal[i] = rotation[i][i_plane]
00525 
00526     return normal
00527 
00528 
00529 def get_bound_project(bound_box, dir):
00530     """Get bounds projection"""
00531     bound_points = [[bound_box[0], bound_box[2], bound_box[4]],
00532                     [bound_box[1], bound_box[2], bound_box[4]],
00533                     [bound_box[0], bound_box[3], bound_box[4]],
00534                     [bound_box[1], bound_box[3], bound_box[4]],
00535                     [bound_box[0], bound_box[2], bound_box[5]],
00536                     [bound_box[1], bound_box[2], bound_box[5]],
00537                     [bound_box[0], bound_box[3], bound_box[5]],
00538                     [bound_box[1], bound_box[3], bound_box[5]]]
00539 
00540     bound_prj = [0, 0, 0]
00541     bound_prj[0] = dot_product(dir, bound_points[0])
00542     bound_prj[1] = bound_prj[0]
00543 
00544     for i in xrange(1, 8):
00545         tmp = dot_product(dir, bound_points[i])
00546         if bound_prj[1] < tmp:
00547             bound_prj[1] = tmp
00548         if bound_prj[0] > tmp:
00549             bound_prj[0] = tmp
00550 
00551     bound_prj[2] = bound_prj[1] - bound_prj[0]
00552     bound_prj[1] = bound_prj[0] + (1.0 - EPS) * bound_prj[2]
00553     bound_prj[0] = bound_prj[0] + EPS * bound_prj[2]
00554     bound_prj[2] = bound_prj[1] - bound_prj[0]
00555 
00556     return bound_prj
00557 
00558 
00559 def get_positions(nb_planes, dir, bounds, displacement):
00560     """Compute plane positions."""
00561     positions = []
00562     bound_prj = get_bound_project(bounds, dir)
00563     if nb_planes > 1:
00564         step = bound_prj[2] / (nb_planes - 1)
00565         abs_displacement = step * displacement
00566         start_pos = bound_prj[0] - 0.5 * step + abs_displacement
00567         for i in xrange(nb_planes):
00568             pos = start_pos + i * step
00569             positions.append(pos)
00570     else:
00571         pos = bound_prj[0] + bound_prj[2] * displacement
00572         positions.append(pos)
00573 
00574     return positions
00575 
00576 
00577 def get_contours(scalar_range, nb_contours):
00578     """Generate contour values."""
00579     contours = []
00580     for i in xrange(nb_contours):
00581         pos = scalar_range[0] + i * (
00582             scalar_range[1] - scalar_range[0]) / (nb_contours - 1)
00583         contours.append(pos)
00584 
00585     return contours
00586 
00587 
00588 def get_nb_components(proxy, entity, field_name):
00589     """Return number of components for the field."""
00590     entity_data_info = None
00591     field_data = proxy.GetFieldDataInformation()
00592 
00593     if field_name in field_data.keys():
00594         entity_data_info = field_data
00595     elif entity == EntityType.CELL:
00596         entity_data_info = proxy.GetCellDataInformation()
00597     elif entity == EntityType.NODE:
00598         entity_data_info = proxy.GetPointDataInformation()
00599 
00600     nb_comp = None
00601     if field_name in entity_data_info.keys():
00602         nb_comp = entity_data_info[field_name].GetNumberOfComponents()
00603     else:
00604         pv_entity = EntityType.get_pvtype(entity)
00605         raise ValueError("Field " + field_name +
00606                          " is unknown for " + pv_entity + "!")
00607 
00608     return nb_comp
00609 
00610 
00611 def get_scale_factor(proxy):
00612     """Compute scale factor."""
00613     if not proxy:
00614         return 0.0
00615 
00616     proxy.UpdatePipeline()
00617     data_info = proxy.GetDataInformation()
00618 
00619     nb_cells = data_info.GetNumberOfCells()
00620     nb_points = data_info.GetNumberOfPoints()
00621     nb_elements = nb_cells if nb_cells > 0  else nb_points
00622     bounds = get_bounds(proxy)
00623 
00624     volume = 1
00625     vol = dim = 0
00626 
00627     for i in xrange(0, 6, 2):
00628         vol = abs(bounds[i + 1] - bounds[i])
00629         if vol > 0:
00630             dim += 1
00631             volume *= vol
00632 
00633     if nb_elements == 0 or dim < 1 / VTK_LARGE_FLOAT:
00634         return 0
00635 
00636     volume /= nb_elements
00637 
00638     return pow(volume, 1 / dim)
00639 
00640 
00641 def get_default_scale(prs_type, proxy, entity, field_name):
00642     """Get default scale factor."""
00643     data_range = get_data_range(proxy, entity, field_name)
00644 
00645     if prs_type == PrsTypeEnum.DEFORMEDSHAPE:
00646         EPS = 1.0 / VTK_LARGE_FLOAT
00647         if abs(data_range[1]) > EPS:
00648             scale_factor = get_scale_factor(proxy)
00649             return scale_factor / data_range[1]
00650     elif prs_type == PrsTypeEnum.PLOT3D:
00651         bounds = get_bounds(proxy)
00652         length = sqrt((bounds[1] - bounds[0]) ** 2 +
00653                       (bounds[3] - bounds[2]) ** 2 +
00654                       (bounds[5] - bounds[4]) ** 2)
00655 
00656         EPS = 0.3
00657         if data_range[1] > 0:
00658             return length / data_range[1] * EPS
00659 
00660     return 0
00661 
00662 
00663 def get_calc_magnitude(proxy, array_entity, array_name):
00664     """Compute magnitude for the given vector array via Calculator.
00665 
00666     Returns:
00667       the calculator object.
00668 
00669     """
00670     calculator = None
00671 
00672     # Transform vector array to scalar array if possible
00673     nb_components = get_nb_components(proxy, array_entity, array_name)
00674     if (nb_components > 1):
00675         calculator = pv.Calculator(proxy)
00676         attribute_mode = "point_data"
00677         if array_entity != EntityType.NODE:
00678             attribute_mode = "cell_data"
00679         calculator.AttributeMode = attribute_mode
00680         if (nb_components == 2):
00681             # Workaroud: calculator unable to compute magnitude
00682             # if number of components equal to 2
00683             func = "sqrt(" + array_name + "_X^2+" + array_name + "_Y^2)"
00684             calculator.Function = func
00685         else:
00686             calculator.Function = "mag(" + array_name + ")"
00687         calculator.ResultArrayName = array_name + "_magnitude"
00688         calculator.UpdatePipeline()
00689 
00690     return calculator
00691 
00692 
00693 def get_add_component_calc(proxy, array_entity, array_name):
00694     """Creates 3-component array from 2-component.
00695 
00696     The first two components is from the original array. The 3rd component
00697     is zero.
00698     If the number of components is not equal to 2 - return original array name.
00699 
00700     Returns:
00701       the calculator object.
00702 
00703     """
00704     calculator = None
00705 
00706     nb_components = get_nb_components(proxy, array_entity, array_name)
00707     if nb_components == 2:
00708         calculator = pv.Calculator(proxy)
00709         attribute_mode = "point_data"
00710         if array_entity != EntityType.NODE:
00711             attribute_mode = "cell_data"
00712         calculator.AttributeMode = attribute_mode
00713         expression = "iHat * " + array_name + "_X + jHat * " + array_name + "_Y + kHat * 0"
00714         calculator.Function = expression
00715         calculator.ResultArrayName = array_name + "_3c"
00716         calculator.UpdatePipeline()
00717 
00718     return calculator
00719 
00720 
00721 def select_all_cells(proxy):
00722     """Select all cell types.
00723 
00724     Used in creation of mesh/submesh presentation.
00725 
00726     """
00727     ### Old API all_cell_types = proxy.CellTypes.Available
00728     all_cell_types = proxy.Entity.Available
00729     ### Old API proxy.CellTypes = all_cell_types
00730     proxy.Entity = all_cell_types
00731     proxy.UpdatePipeline()
00732 
00733 
00734 def select_cells_with_data(proxy, on_points=None, on_cells=None):
00735     """Select cell types with data.
00736 
00737     Only cell types with data for the given fields will be selected.
00738     If no fields defined (neither on points nor on cells) only cell
00739     types with data for even one field (from available) will be selected.
00740 
00741     """
00742     #all_cell_types = proxy.CellTypes.Available
00743     all_cell_types = proxy.Entity.Available
00744     all_arrays = list(proxy.CellArrays.GetData())
00745     all_arrays.extend(proxy.PointArrays.GetData())
00746 
00747     if not all_arrays:
00748         file_name = proxy.FileName.split(os.sep)[-1]
00749         print "Warning: " + file_name + " doesn't contain any data array."
00750 
00751     # List of cell types to be selected
00752     cell_types_on = []
00753 
00754     for cell_type in all_cell_types:
00755         #proxy.CellTypes = [cell_type]
00756         proxy.Entity = [cell_type]
00757         proxy.UpdatePipeline()
00758 
00759         cell_arrays = proxy.GetCellDataInformation().keys()
00760         point_arrays = proxy.GetPointDataInformation().keys()
00761 
00762         if on_points or on_cells:
00763             if on_points is None:
00764                 on_points = []
00765             if on_cells is None:
00766                 on_cells = []
00767 
00768             if (all(array in cell_arrays for array in on_cells) and
00769                 all(array in point_arrays for array in on_points)):
00770                 # Add cell type to the list
00771                 cell_types_on.append(cell_type)
00772         else:
00773             in_arrays = lambda array: ((array in cell_arrays) or
00774                                        (array in point_arrays))
00775             if any(in_arrays(array) for array in all_arrays):
00776                 cell_types_on.append(cell_type)
00777 
00778     # Select cell types
00779     #proxy.CellTypes = cell_types_on
00780     proxy.Entity = cell_types_on
00781     proxy.UpdatePipeline()
00782 
00783 
00784 def extract_groups_for_field(proxy, field_name, field_entity, force=False):
00785     """Exctract only groups which have the field.
00786 
00787     Arguments:
00788       proxy: the pipeline object, containig data
00789       field_name: the field name
00790       field_entity: the field entity
00791       force: if True - ExtractGroup object will be created in any case
00792 
00793     Returns:
00794       ExtractGroup object: if not all groups have the field or
00795       the force argument is true
00796       The initial proxy: if no groups had been filtered.
00797 
00798     """
00799     source = proxy
00800 
00801     # Remember the state
00802     initial_groups = list(proxy.Groups)
00803 
00804     # Get data information for the field entity
00805     entity_data_info = None
00806     field_data = proxy.GetFieldDataInformation()
00807 
00808     if field_name in field_data.keys():
00809         entity_data_info = field_data
00810     elif field_entity == EntityType.CELL:
00811         entity_data_info = proxy.GetCellDataInformation()
00812     elif field_entity == EntityType.NODE:
00813         entity_data_info = proxy.GetPointDataInformation()
00814 
00815     # Collect groups for extraction
00816     groups_to_extract = []
00817 
00818     for group in initial_groups:
00819         proxy.Groups = [group]
00820         proxy.UpdatePipeline()
00821         if field_name in entity_data_info.keys():
00822             groups_to_extract.append(group)
00823 
00824     # Restore state
00825     proxy.Groups = initial_groups
00826     proxy.UpdatePipeline()
00827 
00828     # Extract groups if necessary
00829     if force or (len(groups_to_extract) < len(initial_groups)):
00830         extract_group = pv.ExtractGroup(proxy)
00831         extract_group.Groups = groups_to_extract
00832         extract_group.UpdatePipeline()
00833         source = extract_group
00834 
00835     return source
00836 
00837 
00838 def if_possible(proxy, field_name, entity, prs_type):
00839     """Check if the presentation creation is possible on the given field."""
00840     result = True
00841     if (prs_type == PrsTypeEnum.DEFORMEDSHAPE or
00842         prs_type == PrsTypeEnum.DEFORMEDSHAPESCALARMAP or
00843         prs_type == PrsTypeEnum.VECTORS or
00844         prs_type == PrsTypeEnum.STREAMLINES):
00845         nb_comp = get_nb_components(proxy, entity, field_name)
00846         result = (nb_comp > 1)
00847     elif (prs_type == PrsTypeEnum.GAUSSPOINTS):
00848         result = (entity == EntityType.CELL or
00849                   field_name in proxy.QuadraturePointArrays.Available)
00850     elif (prs_type == PrsTypeEnum.MESH):
00851         result = len(get_group_names(proxy, field_name, entity)) > 0
00852 
00853     return result
00854 
00855 
00856 def add_scalar_bar(field_name, nb_components,
00857                    vector_mode, lookup_table, time_value):
00858     """Add scalar bar with predefined properties."""
00859     global _current_bar
00860 
00861     # Construct bar title
00862     title = "\n".join([field_name, str(time_value)])
00863     if nb_components > 1:
00864         title = "\n".join([title, vector_mode])
00865 
00866     # Create scalar bar
00867     scalar_bar = pv.CreateScalarBar(Enabled=1)
00868     scalar_bar.Orientation = 'Vertical'
00869     scalar_bar.Title = title
00870     scalar_bar.LookupTable = lookup_table
00871 
00872     # Set default properties same as in Post-Pro
00873     scalar_bar.NumberOfLabels = 5
00874     scalar_bar.AutomaticLabelFormat = 0
00875     scalar_bar.LabelFormat = '%-#6.6g'
00876     # Title
00877     scalar_bar.TitleFontFamily = 'Arial'
00878     scalar_bar.TitleFontSize = 8
00879     scalar_bar.TitleBold = 1
00880     scalar_bar.TitleItalic = 1
00881     scalar_bar.TitleShadow = 1
00882     # Labels
00883     scalar_bar.LabelFontFamily = 'Arial'
00884     scalar_bar.LabelFontSize = 8
00885     scalar_bar.LabelBold = 1
00886     scalar_bar.LabelItalic = 1
00887     scalar_bar.LabelShadow = 1
00888 
00889     # Add the scalar bar to the view
00890     pv.GetRenderView().Representations.append(scalar_bar)
00891 
00892     # Reassign the current bar
00893     _current_bar = scalar_bar
00894 
00895     return scalar_bar
00896 
00897 
00898 def get_bar():
00899     """Get current scalar bar."""
00900     global _current_bar
00901 
00902     return _current_bar
00903 
00904 
00905 def get_lookup_table(field_name, nb_components, vector_mode='Magnitude'):
00906     """Get lookup table for the given field."""
00907     lookup_table = pv.GetLookupTableForArray(field_name, nb_components)
00908 
00909     if vector_mode == 'Magnitude':
00910         lookup_table.VectorMode = vector_mode
00911     elif vector_mode == 'X':
00912         lookup_table.VectorMode = 'Component'
00913         lookup_table.VectorComponent = 0
00914     elif vector_mode == 'Y':
00915         lookup_table.VectorMode = 'Component'
00916         lookup_table.VectorComponent = 1
00917     elif vector_mode == 'Z':
00918         lookup_table.VectorMode = 'Component'
00919         lookup_table.VectorComponent = 2
00920     else:
00921         raise ValueError("Incorrect vector mode: " + vector_mode)
00922 
00923     lookup_table.Discretize = 0
00924     lookup_table.ColorSpace = 'HSV'
00925     lookup_table.LockScalarRange = 0
00926 
00927     return lookup_table
00928 
00929 
00930 def get_group_mesh_name(full_group_name):
00931     """Return mesh name of the group by its full name."""
00932     aList = full_group_name.split('/')
00933     if len(aList) >= 2 :
00934         group_name = full_group_name.split('/')[1]
00935         return group_name
00936 
00937 
00938 def get_group_entity(full_group_name):
00939     """Return entity type of the group by its full name."""
00940     aList = full_group_name.split('/')
00941     if len(aList) >= 3 :
00942         entity_name = full_group_name.split('/')[2]
00943         entity = EntityType.get_type(entity_name)
00944         return entity
00945 
00946 
00947 def get_group_short_name(full_group_name):
00948     """Return short name of the group by its full name."""
00949     aList = full_group_name.split('/')
00950     if len(aList) >= 4 :
00951         short_name = full_group_name.split('/')[3]
00952         return short_name
00953 
00954 
00955 def get_mesh_names(proxy):
00956     """Return all mesh names in the given proxy as a set."""
00957     groups = proxy.Groups.Available
00958     mesh_names = set([get_group_mesh_name(item) for item in groups])
00959 
00960     return mesh_names
00961 
00962 
00963 def get_group_names(proxy, mesh_name, entity, wo_nogroups=False):
00964     """Return full names of all groups of the given entity type
00965     from the mesh with the given name as a list.
00966     """
00967     groups = proxy.Groups.Available
00968 
00969     condition = lambda item: (get_group_mesh_name(item) == mesh_name and
00970                               get_group_entity(item) == entity)
00971     group_names = [item for item in groups if condition(item)]
00972 
00973     if wo_nogroups:
00974         # Remove "No_Group" group
00975         not_no_group = lambda item: get_group_short_name(item) != "No_Group"
00976         group_names = filter(not_no_group, group_names)
00977 
00978     return group_names
00979 
00980 
00981 def get_time(proxy, timestamp_nb):
00982     """Get time value by timestamp number."""
00983     # Check timestamp number
00984     timestamps = proxy.TimestepValues.GetData()
00985     if ((timestamp_nb - 1) not in xrange(len(timestamps))):
00986         raise ValueError("Timestamp number is out of range: " + timestamp_nb)
00987 
00988     # Return time value
00989     return timestamps[timestamp_nb - 1]
00990 
00991 
00992 def create_prs(prs_type, proxy, field_entity, field_name, timestamp_nb):
00993     """Auxiliary function.
00994 
00995     Build presentation of the given type on the given field and
00996     timestamp number.
00997     Set the presentation properties like visu.CreatePrsForResult() do.
00998 
00999     """
01000     prs = None
01001 
01002     if prs_type == PrsTypeEnum.SCALARMAP:
01003         prs = ScalarMapOnField(proxy, field_entity, field_name, timestamp_nb)
01004     elif prs_type == PrsTypeEnum.CUTPLANES:
01005         prs = CutPlanesOnField(proxy, field_entity, field_name, timestamp_nb,
01006                                orientation=Orientation.ZX)
01007     elif prs_type == PrsTypeEnum.CUTLINES:
01008         prs = CutLinesOnField(proxy, field_entity, field_name, timestamp_nb,
01009                               orientation1=Orientation.XY,
01010                               orientation2=Orientation.ZX)
01011     elif prs_type == PrsTypeEnum.DEFORMEDSHAPE:
01012         prs = DeformedShapeOnField(proxy, field_entity,
01013                                    field_name, timestamp_nb)
01014     elif prs_type == PrsTypeEnum.DEFORMEDSHAPESCALARMAP:
01015         prs = DeformedShapeAndScalarMapOnField(proxy, field_entity,
01016                                                field_name, timestamp_nb)
01017     elif prs_type == PrsTypeEnum.VECTORS:
01018         prs = VectorsOnField(proxy, field_entity, field_name, timestamp_nb)
01019     elif prs_type == PrsTypeEnum.PLOT3D:
01020         prs = Plot3DOnField(proxy, field_entity, field_name, timestamp_nb)
01021     elif prs_type == PrsTypeEnum.ISOSURFACES:
01022         prs = IsoSurfacesOnField(proxy, field_entity, field_name, timestamp_nb)
01023     elif prs_type == PrsTypeEnum.GAUSSPOINTS:
01024         prs = GaussPointsOnField(proxy, field_entity, field_name, timestamp_nb)
01025     elif prs_type == PrsTypeEnum.STREAMLINES:
01026         prs = StreamLinesOnField(proxy, field_entity, field_name, timestamp_nb)
01027     else:
01028         raise ValueError("Unexistent presentation type.")
01029 
01030     return prs
01031 
01032 
01033 # Functions for building Post-Pro presentations
01034 def ScalarMapOnField(proxy, entity, field_name, timestamp_nb,
01035                      vector_mode='Magnitude'):
01036     """Creates Scalar Map presentation on the given field.
01037 
01038     Arguments:
01039       proxy: the pipeline object, containig data
01040       entity: the entity type from PrsTypeEnum
01041       field_name: the field name
01042       timestamp_nb: the number of time step (1, 2, ...)
01043       vector_mode: the mode of transformation of vector values
01044       into scalar values, applicable only if the field contains vector values.
01045       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01046 
01047     Returns:
01048       Scalar Map as representation object.
01049 
01050     """
01051     # We don't need mesh parts with no data on them
01052     if entity == EntityType.NODE:
01053         select_cells_with_data(proxy, on_points=[field_name])
01054     else:
01055         select_cells_with_data(proxy, on_cells=[field_name])
01056 
01057     # Check vector mode
01058     nb_components = get_nb_components(proxy, entity, field_name)
01059     check_vector_mode(vector_mode, nb_components)
01060 
01061     # Get time value
01062     time_value = get_time(proxy, timestamp_nb)
01063 
01064     # Set timestamp
01065     pv.GetRenderView().ViewTime = time_value
01066     pv.UpdatePipeline(time_value, proxy)
01067 
01068     # Extract only groups with data for the field
01069     new_proxy = extract_groups_for_field(proxy, field_name, entity,
01070                                          force=True)
01071 
01072     # Get Scalar Map representation object
01073     scalarmap = pv.GetRepresentation(new_proxy)
01074 
01075     # Get lookup table
01076     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01077 
01078     # Set field range if necessary
01079     data_range = get_data_range(proxy, entity,
01080                                 field_name, vector_mode)
01081     lookup_table.LockScalarRange = 1
01082     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01083     # Set properties
01084     scalarmap.ColorAttributeType = EntityType.get_pvtype(entity)
01085     scalarmap.ColorArrayName = field_name
01086     scalarmap.LookupTable = lookup_table
01087 
01088     # Add scalar bar
01089     bar_title = field_name + ", " + str(time_value)
01090     if (nb_components > 1):
01091         bar_title += "\n" + vector_mode
01092     add_scalar_bar(field_name, nb_components, vector_mode,
01093                    lookup_table, time_value)
01094 
01095     return scalarmap
01096 
01097 
01098 def CutPlanesOnField(proxy, entity, field_name, timestamp_nb,
01099                      nb_planes=10, orientation=Orientation.YZ,
01100                      angle1=0, angle2=0,
01101                      displacement=0.5, vector_mode='Magnitude'):
01102     """Creates Cut Planes presentation on the given field.
01103 
01104     Arguments:
01105       proxy: the pipeline object, containig data
01106       entity: the entity type from PrsTypeEnum
01107       field_name: the field name
01108       timestamp_nb: the number of time step (1, 2, ...)
01109       nb_planes: number of cutting planes
01110       orientation: cutting planes orientation in 3D space
01111       angle1: rotation of the planes in 3d space around the first axis of the
01112       selected orientation (X axis for XY, Y axis for YZ, Z axis for ZX).
01113       The angle of rotation is set in degrees. Acceptable range: [-45, 45].
01114       angle2: rotation of the planes in 3d space around the second axis of the
01115       selected orientation. Acceptable range: [-45, 45].
01116       displacement: the displacement of the planes into one or another side
01117       vector_mode: the mode of transformation of vector values
01118       into scalar values, applicable only if the field contains vector values.
01119       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01120 
01121     Returns:
01122       Cut Planes as representation object.
01123 
01124     """
01125     # Check vector mode
01126     nb_components = get_nb_components(proxy, entity, field_name)
01127     check_vector_mode(vector_mode, nb_components)
01128 
01129     # Get time value
01130     time_value = get_time(proxy, timestamp_nb)
01131 
01132     # Set timestamp
01133     pv.GetRenderView().ViewTime = time_value
01134     pv.UpdatePipeline(time_value, proxy)
01135 
01136     # Create slice filter
01137     slice_filter = pv.Slice(proxy)
01138     slice_filter.SliceType = "Plane"
01139 
01140     # Set cut planes normal
01141     normal = get_normal_by_orientation(orientation,
01142                                        radians(angle1), radians(angle2))
01143     slice_filter.SliceType.Normal = normal
01144 
01145     # Set cut planes positions
01146     positions = get_positions(nb_planes, normal,
01147                               get_bounds(proxy), displacement)
01148     slice_filter.SliceOffsetValues = positions
01149 
01150     # Get Cut Planes representation object
01151     cut_planes = pv.GetRepresentation(slice_filter)
01152 
01153     # Get lookup table
01154     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01155 
01156     # Set field range if necessary
01157     data_range = get_data_range(proxy, entity,
01158                                 field_name, vector_mode)
01159     lookup_table.LockScalarRange = 1
01160     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01161 
01162     # Set properties
01163     cut_planes.ColorAttributeType = EntityType.get_pvtype(entity)
01164     cut_planes.ColorArrayName = field_name
01165     cut_planes.LookupTable = lookup_table
01166 
01167     # Add scalar bar
01168     add_scalar_bar(field_name, nb_components,
01169                    vector_mode, lookup_table, time_value)
01170 
01171     return cut_planes
01172 
01173 
01174 def CutLinesOnField(proxy, entity, field_name, timestamp_nb,
01175                     nb_lines=10,
01176                     orientation1=Orientation.XY,
01177                     base_angle1=0, base_angle2=0,
01178                     orientation2=Orientation.YZ,
01179                     cut_angle1=0, cut_angle2=0,
01180                     displacement1=0.5, displacement2=0.5,
01181                     generate_curves=False,
01182                     vector_mode='Magnitude'):
01183     """Creates Cut Lines presentation on the given field.
01184 
01185     Arguments:
01186       proxy: the pipeline object, containig data
01187       entity: the entity type from PrsTypeEnum
01188       field_name: the field name
01189       timestamp_nb: the number of time step (1, 2, ...)
01190       nb_lines: number of lines
01191       orientation1: base plane orientation in 3D space
01192       base_angle1: rotation of the base plane in 3d space around the first
01193       axis of the orientation1 (X axis for XY, Y axis for YZ, Z axis for ZX).
01194       The angle of rotation is set in degrees. Acceptable range: [-45, 45].
01195       base_angle2: rotation of the base plane in 3d space around the second
01196       axis of the orientation1. Acceptable range: [-45, 45].
01197       orientation2: cutting planes orientation in 3D space
01198       cut_angle1: rotation of the cut planes in 3d space around the first
01199       axis of the orientation2. Acceptable range: [-45, 45].
01200       cut_angle2: rotation of the cuting planes in 3d space around the second
01201       axis of the orientation2. Acceptable range: [-45, 45].
01202       displacement1: base plane displacement
01203       displacement2: cutting planes displacement
01204       generate_curves: if true, 'PlotOverLine' filter will be created
01205       for each cut line
01206       vector_mode: the mode of transformation of vector values
01207       into scalar values, applicable only if the field contains vector values.
01208       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01209 
01210     Returns:
01211       Cut Lines as representation object if generate_curves == False,
01212       (Cut Lines as representation object, list of 'PlotOverLine') otherwise
01213 
01214     """
01215     # Check vector mode
01216     nb_components = get_nb_components(proxy, entity, field_name)
01217     check_vector_mode(vector_mode, nb_components)
01218 
01219     # Get time value
01220     time_value = get_time(proxy, timestamp_nb)
01221 
01222     # Set timestamp
01223     pv.GetRenderView().ViewTime = time_value
01224     pv.UpdatePipeline(time_value, proxy)
01225 
01226     # Create base plane
01227     base_plane = pv.Slice(proxy)
01228     base_plane.SliceType = "Plane"
01229 
01230     # Set base plane normal
01231     base_normal = get_normal_by_orientation(orientation1,
01232                                             radians(base_angle1),
01233                                             radians(base_angle2))
01234     base_plane.SliceType.Normal = base_normal
01235 
01236     # Set base plane position
01237     base_position = get_positions(1, base_normal,
01238                                   get_bounds(proxy), displacement1)
01239     base_plane.SliceOffsetValues = base_position
01240 
01241     # Check base plane
01242     base_plane.UpdatePipeline()
01243     if (base_plane.GetDataInformation().GetNumberOfCells() == 0):
01244         base_plane = proxy
01245 
01246     # Create cutting planes
01247     cut_planes = pv.Slice(base_plane)
01248     cut_planes.SliceType = "Plane"
01249 
01250     # Set cutting planes normal and get positions
01251     cut_normal = get_normal_by_orientation(orientation2,
01252                                            radians(cut_angle1),
01253                                            radians(cut_angle2))
01254     cut_planes.SliceType.Normal = cut_normal
01255 
01256     # Set cutting planes position
01257     cut_positions = get_positions(nb_lines, cut_normal,
01258                                   get_bounds(base_plane), displacement2)
01259 
01260     # Generate curves
01261     curves = []
01262     if generate_curves:
01263         index = 0
01264         for pos in cut_positions:
01265             # Get points for plot over line objects
01266             cut_planes.SliceOffsetValues = pos
01267             cut_planes.UpdatePipeline()
01268             bounds = get_bounds(cut_planes)
01269             point1 = [bounds[0], bounds[2], bounds[4]]
01270             point2 = [bounds[1], bounds[3], bounds[5]]
01271 
01272             # Create plot over line filter
01273             pol = pv.PlotOverLine(cut_planes,
01274                                   Source="High Resolution Line Source")
01275             pv.RenameSource('Y' + str(index), pol)
01276             pol.Source.Point1 = point1
01277             pol.Source.Point2 = point2
01278             pol.UpdatePipeline()
01279             curves.append(pol)
01280 
01281             index += 1
01282 
01283     cut_planes.SliceOffsetValues = cut_positions
01284     cut_planes.UpdatePipeline()
01285 
01286     # Get Cut Lines representation object
01287     cut_lines = pv.GetRepresentation(cut_planes)
01288 
01289     # Get lookup table
01290     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01291 
01292     # Set field range if necessary
01293     data_range = get_data_range(proxy, entity,
01294                                 field_name, vector_mode)
01295     lookup_table.LockScalarRange = 1
01296     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01297 
01298     # Set properties
01299     cut_lines.ColorAttributeType = EntityType.get_pvtype(entity)
01300     cut_lines.ColorArrayName = field_name
01301     cut_lines.LookupTable = lookup_table
01302 
01303     # Set wireframe represenatation mode
01304     cut_lines.Representation = 'Wireframe'
01305 
01306     # Add scalar bar
01307     add_scalar_bar(field_name, nb_components,
01308                    vector_mode, lookup_table, time_value)
01309 
01310     result = cut_lines
01311     # If curves were generated return tuple (cut lines, list of curves)
01312     if curves:
01313         result = cut_lines, curves
01314 
01315     return result
01316 
01317 
01318 def CutSegmentOnField(proxy, entity, field_name, timestamp_nb,
01319                       point1, point2, vector_mode='Magnitude'):
01320     """Creates Cut Segment presentation on the given field.
01321 
01322     Arguments:
01323       proxy: the pipeline object, containig data
01324       entity: the entity type from PrsTypeEnum
01325       field_name: the field name
01326       timestamp_nb: the number of time step (1, 2, ...)
01327       point1: set the first point of the segment (as [x, y, z])
01328       point1: set the second point of the segment (as [x, y, z])
01329       vector_mode: the mode of transformation of vector values
01330       into scalar values, applicable only if the field contains vector values.
01331       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01332 
01333     Returns:
01334       Cut Segment as 3D representation object.
01335 
01336     """
01337     # Check vector mode
01338     nb_components = get_nb_components(proxy, entity, field_name)
01339     check_vector_mode(vector_mode, nb_components)
01340 
01341     # Get time value
01342     time_value = get_time(proxy, timestamp_nb)
01343 
01344     # Set timestamp
01345     pv.GetRenderView().ViewTime = time_value
01346     pv.UpdatePipeline(time_value, proxy)
01347 
01348     # Create plot over line filter
01349     pol = pv.PlotOverLine(proxy, Source="High Resolution Line Source")
01350     pol.Source.Point1 = point1
01351     pol.Source.Point2 = point2
01352     pol.UpdatePipeline()
01353 
01354     # Get Cut Segment representation object
01355     cut_segment = pv.GetRepresentation(pol)
01356 
01357     # Get lookup table
01358     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01359 
01360     # Set field range if necessary
01361     data_range = get_data_range(proxy, entity,
01362                                 field_name, vector_mode)
01363     lookup_table.LockScalarRange = 1
01364     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01365 
01366     # Set properties
01367     cut_segment.ColorAttributeType = EntityType.get_pvtype(entity)
01368     cut_segment.ColorArrayName = field_name
01369     cut_segment.LookupTable = lookup_table
01370 
01371     # Set wireframe represenatation mode
01372     cut_segment.Representation = 'Wireframe'
01373 
01374     # Add scalar bar
01375     add_scalar_bar(field_name, nb_components,
01376                    vector_mode, lookup_table, time_value)
01377 
01378     return cut_segment
01379 
01380 
01381 def VectorsOnField(proxy, entity, field_name, timestamp_nb,
01382                    scale_factor=None,
01383                    glyph_pos=GlyphPos.TAIL, glyph_type='2D Glyph',
01384                    is_colored=False, vector_mode='Magnitude'):
01385     """Creates Vectors presentation on the given field.
01386 
01387     Arguments:
01388       proxy: the pipeline object, containig data
01389       entity: the entity type from PrsTypeEnum
01390       field_name: the field name
01391       timestamp_nb: the number of time step (1, 2, ...)
01392       scale_factor: scale factor
01393       glyph_pos: the position of glyphs
01394       glyph_type: the type of glyphs
01395       is_colored: this option allows to color the presentation according to
01396       the corresponding data array values
01397       vector_mode: the mode of transformation of vector values
01398       into scalar values, applicable only if the field contains vector values.
01399       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01400 
01401     Returns:
01402       Vectors as representation object.
01403 
01404     """
01405     # Check vector mode
01406     nb_components = get_nb_components(proxy, entity, field_name)
01407     check_vector_mode(vector_mode, nb_components)
01408 
01409     # Get time value
01410     time_value = get_time(proxy, timestamp_nb)
01411 
01412     # Set timestamp
01413     pv.GetRenderView().ViewTime = time_value
01414     pv.UpdatePipeline(time_value, proxy)
01415 
01416     # Extract only groups with data for the field
01417     new_proxy = extract_groups_for_field(proxy, field_name, entity)
01418     source = new_proxy
01419 
01420     # Cell centers
01421     if is_data_on_cells(proxy, field_name):
01422         cell_centers = pv.CellCenters(source)
01423         cell_centers.VertexCells = 1
01424         source = cell_centers
01425 
01426     vector_array = field_name
01427     # If the given vector array has only 2 components, add the third one
01428     if nb_components == 2:
01429         calc = get_add_component_calc(source, EntityType.NODE, field_name)
01430         vector_array = calc.ResultArrayName
01431         source = calc
01432 
01433     # Glyph
01434     glyph = pv.Glyph(source)
01435     glyph.Vectors = vector_array
01436     glyph.ScaleMode = 'vector'
01437     glyph.MaskPoints = 0
01438 
01439     # Set glyph type
01440     glyph.GlyphType = glyph_type
01441     if glyph_type == '2D Glyph':
01442         glyph.GlyphType.GlyphType = 'Arrow'
01443     elif glyph_type == 'Cone':
01444         glyph.GlyphType.Resolution = 7
01445         glyph.GlyphType.Height = 2
01446         glyph.GlyphType.Radius = 0.2
01447 
01448     # Set glyph position if possible
01449     if glyph.GlyphType.GetProperty("Center"):
01450         if (glyph_pos == GlyphPos.TAIL):
01451             glyph.GlyphType.Center = [0.5, 0.0, 0.0]
01452         elif (glyph_pos == GlyphPos.HEAD):
01453             glyph.GlyphType.Center = [-0.5, 0.0, 0.0]
01454         elif (glyph_pos == GlyphPos.CENTER):
01455             glyph.GlyphType.Center = [0.0, 0.0, 0.0]
01456 
01457     if scale_factor is not None:
01458         glyph.SetScaleFactor = scale_factor
01459     else:
01460         def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE,
01461                                       new_proxy, entity, field_name)
01462         glyph.SetScaleFactor = def_scale
01463 
01464     glyph.UpdatePipeline()
01465 
01466     # Get Vectors representation object
01467     vectors = pv.GetRepresentation(glyph)
01468 
01469     # Get lookup table
01470     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01471 
01472     # Set field range if necessary
01473     data_range = get_data_range(proxy, entity,
01474                                 field_name, vector_mode)
01475     lookup_table.LockScalarRange = 1
01476     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01477 
01478     # Set properties
01479     if (is_colored):
01480         vectors.ColorArrayName = 'GlyphVector'
01481     else:
01482         vectors.ColorArrayName = ''
01483     vectors.LookupTable = lookup_table
01484 
01485     vectors.LineWidth = 1.0
01486 
01487     # Set wireframe represenatation mode
01488     vectors.Representation = 'Wireframe'
01489 
01490     # Add scalar bar
01491     add_scalar_bar(field_name, nb_components,
01492                    vector_mode, lookup_table, time_value)
01493 
01494     return vectors
01495 
01496 
01497 def DeformedShapeOnField(proxy, entity, field_name,
01498                          timestamp_nb,
01499                          scale_factor=None, is_colored=False,
01500                          vector_mode='Magnitude'):
01501     """Creates Defromed Shape presentation on the given field.
01502 
01503     Arguments:
01504       proxy: the pipeline object, containig data
01505       entity: the entity type from PrsTypeEnum
01506       field_name: the field name
01507       timestamp_nb: the number of time step (1, 2, ...)
01508       scale_factor: scale factor of the deformation
01509       is_colored: this option allows to color the presentation according to
01510       the corresponding data array values
01511       vector_mode: the mode of transformation of vector values
01512       into scalar values, applicable only if the field contains vector values.
01513       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01514 
01515     Returns:
01516       Defromed Shape as representation object.
01517 
01518     """
01519     # We don't need mesh parts with no data on them
01520     if entity == EntityType.NODE:
01521         select_cells_with_data(proxy, on_points=[field_name])
01522     else:
01523         select_cells_with_data(proxy, on_cells=[field_name])
01524 
01525     # Check vector mode
01526     nb_components = get_nb_components(proxy, entity, field_name)
01527     check_vector_mode(vector_mode, nb_components)
01528 
01529     # Get time value
01530     time_value = get_time(proxy, timestamp_nb)
01531 
01532     # Set timestamp
01533     pv.GetRenderView().ViewTime = time_value
01534     pv.UpdatePipeline(time_value, proxy)
01535 
01536     # Extract only groups with data for the field
01537     new_proxy = extract_groups_for_field(proxy, field_name, entity)
01538 
01539     # Do merge
01540     source = pv.MergeBlocks(new_proxy)
01541 
01542     # Cell data to point data
01543     if is_data_on_cells(proxy, field_name):
01544         cell_to_point = pv.CellDatatoPointData()
01545         cell_to_point.PassCellData = 1
01546         source = cell_to_point
01547 
01548     vector_array = field_name
01549     # If the given vector array has only 2 components, add the third one
01550     if nb_components == 2:
01551         calc = get_add_component_calc(source, EntityType.NODE, field_name)
01552         vector_array = calc.ResultArrayName
01553         source = calc
01554 
01555     # Warp by vector
01556     warp_vector = pv.WarpByVector(source)
01557     warp_vector.Vectors = [vector_array]
01558     if scale_factor is not None:
01559         warp_vector.ScaleFactor = scale_factor
01560     else:
01561         def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE,
01562                                       proxy, entity, field_name)
01563         warp_vector.ScaleFactor = def_scale
01564 
01565     # Get Deformed Shape representation object
01566     defshape = pv.GetRepresentation(warp_vector)
01567 
01568     # Get lookup table
01569     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01570 
01571     # Set field range if necessary
01572     data_range = get_data_range(proxy, entity,
01573                                 field_name, vector_mode)
01574     lookup_table.LockScalarRange = 1
01575     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01576 
01577     # Set properties
01578     if is_colored:
01579         defshape.ColorAttributeType = EntityType.get_pvtype(entity)
01580         defshape.ColorArrayName = field_name
01581     else:
01582         defshape.ColorArrayName = ''
01583     defshape.LookupTable = lookup_table
01584 
01585     # Set wireframe represenatation mode
01586     defshape.Representation = 'Wireframe'
01587 
01588     # Add scalar bar
01589     add_scalar_bar(field_name, nb_components,
01590                    vector_mode, lookup_table, time_value)
01591 
01592     return defshape
01593 
01594 
01595 def DeformedShapeAndScalarMapOnField(proxy, entity, field_name,
01596                                      timestamp_nb,
01597                                      scale_factor=None,
01598                                      scalar_entity=None,
01599                                      scalar_field_name=None,
01600                                      vector_mode='Magnitude'):
01601     """Creates Defromed Shape And Scalar Map presentation on the given field.
01602 
01603     Arguments:
01604       proxy: the pipeline object, containig data
01605       entity: the entity type from PrsTypeEnum
01606       field_name: the field name
01607       timestamp_nb: the number of time step (1, 2, ...)
01608       scale_factor: scale factor of the deformation
01609       scalar_entity: scalar field entity
01610       scalar_field_name: scalar field, i.e. the field for coloring
01611       vector_mode: the mode of transformation of vector values
01612       into scalar values, applicable only if the field contains vector values.
01613       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01614 
01615     Returns:
01616       Defromed Shape And Scalar Map as representation object.
01617 
01618     """
01619     # We don't need mesh parts with no data on them
01620     on_points = []
01621     on_cells = []
01622 
01623     if entity == EntityType.NODE:
01624         on_points.append(field_name)
01625     else:
01626         on_cells.append(field_name)
01627 
01628     if scalar_entity and scalar_field_name:
01629         if scalar_entity == EntityType.NODE:
01630             on_points.append(scalar_field_name)
01631         else:
01632             on_cells.append(scalar_field_name)
01633 
01634     select_cells_with_data(proxy, on_points, on_cells)
01635 
01636     # Check vector mode
01637     nb_components = get_nb_components(proxy, entity, field_name)
01638     check_vector_mode(vector_mode, nb_components)
01639 
01640     # Get time value
01641     time_value = get_time(proxy, timestamp_nb)
01642 
01643     # Set timestamp
01644     pv.GetRenderView().ViewTime = time_value
01645     pv.UpdatePipeline(time_value, proxy)
01646 
01647     # Set scalar field by default
01648     scalar_field_entity = scalar_entity
01649     scalar_field = scalar_field_name
01650     if (scalar_field_entity is None) or (scalar_field is None):
01651         scalar_field_entity = entity
01652         scalar_field = field_name
01653 
01654     # Extract only groups with data for the field
01655     new_proxy = extract_groups_for_field(proxy, field_name, entity)
01656 
01657     # Do merge
01658     source = pv.MergeBlocks(new_proxy)
01659 
01660     # Cell data to point data
01661     if is_data_on_cells(proxy, field_name):
01662         cell_to_point = pv.CellDatatoPointData(source)
01663         cell_to_point.PassCellData = 1
01664         source = cell_to_point
01665 
01666     vector_array = field_name
01667     # If the given vector array has only 2 components, add the third one
01668     if nb_components == 2:
01669         calc = get_add_component_calc(source, EntityType.NODE, field_name)
01670         vector_array = calc.ResultArrayName
01671         source = calc
01672 
01673     # Warp by vector
01674     warp_vector = pv.WarpByVector(source)
01675     warp_vector.Vectors = [vector_array]
01676     if scale_factor is not None:
01677         warp_vector.ScaleFactor = scale_factor
01678     else:
01679         def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE,
01680                                       new_proxy, entity, field_name)
01681         warp_vector.ScaleFactor = def_scale
01682 
01683     # Get Defromed Shape And Scalar Map representation object
01684     defshapemap = pv.GetRepresentation(warp_vector)
01685 
01686     # Get lookup table
01687     lookup_table = get_lookup_table(scalar_field, nb_components, vector_mode)
01688 
01689     # Set field range if necessary
01690     data_range = get_data_range(proxy, scalar_field_entity,
01691                                 scalar_field, vector_mode)
01692     lookup_table.LockScalarRange = 1
01693     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01694 
01695     # Set properties
01696     defshapemap.ColorArrayName = scalar_field
01697     defshapemap.LookupTable = lookup_table
01698     defshapemap.ColorAttributeType = EntityType.get_pvtype(scalar_field_entity)
01699 
01700     # Add scalar bar
01701     add_scalar_bar(field_name, nb_components,
01702                    vector_mode, lookup_table, time_value)
01703 
01704     return defshapemap
01705 
01706 
01707 def Plot3DOnField(proxy, entity, field_name, timestamp_nb,
01708                   orientation=Orientation.AUTO,
01709                   angle1=0, angle2=0,
01710                   position=0.5, is_relative=True,
01711                   scale_factor=None,
01712                   is_contour=False, nb_contours=32,
01713                   vector_mode='Magnitude'):
01714     """Creates Plot 3D presentation on the given field.
01715 
01716     Arguments:
01717       proxy: the pipeline object, containig data
01718       entity: the entity type from PrsTypeEnum
01719       field_name: the field name
01720       timestamp_nb: the number of time step (1, 2, ...)
01721       orientation: the cut plane plane orientation in 3D space, if
01722       the input is planar - will not be taken into account
01723       angle1: rotation of the cut plane in 3d space around the first axis
01724       of the selected orientation (X axis for XY, Y axis for YZ,
01725       Z axis for ZX).
01726       The angle of rotation is set in degrees. Acceptable range: [-45, 45].
01727       angle2: rotation of the cut plane in 3d space around the second axis
01728       of the selected orientation. Acceptable range: [-45, 45].
01729       position: position of the cut plane in the object (ranging from 0 to 1).
01730       The value 0.5 corresponds to cutting by halves.
01731       is_relative: defines if the cut plane position is relative or absolute
01732       scale_factor: deformation scale factor
01733       is_contour: if True - Plot 3D will be represented with a set of contours,
01734       otherwise - Plot 3D will be represented with a smooth surface
01735       nb_contours: number of contours, applied if is_contour is True
01736       vector_mode: the mode of transformation of vector values
01737       into scalar values, applicable only if the field contains vector values.
01738       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01739 
01740     Returns:
01741       Plot 3D as representation object.
01742 
01743     """
01744     # We don't need mesh parts with no data on them
01745     if entity == EntityType.NODE:
01746         select_cells_with_data(proxy, on_points=[field_name])
01747     else:
01748         select_cells_with_data(proxy, on_cells=[field_name])
01749 
01750     # Check vector mode
01751     nb_components = get_nb_components(proxy, entity, field_name)
01752     check_vector_mode(vector_mode, nb_components)
01753 
01754     # Get time value
01755     time_value = get_time(proxy, timestamp_nb)
01756 
01757     # Set timestamp
01758     pv.GetRenderView().ViewTime = time_value
01759     pv.UpdatePipeline(time_value, proxy)
01760 
01761     # Extract only groups with data for the field
01762     new_proxy = extract_groups_for_field(proxy, field_name, entity)
01763 
01764     # Do merge
01765     merge_blocks = pv.MergeBlocks(new_proxy)
01766     merge_blocks.UpdatePipeline()
01767 
01768     poly_data = None
01769 
01770     # Cutting plane
01771 
01772     # Define orientation if necessary (auto mode)
01773     plane_orientation = orientation
01774     if (orientation == Orientation.AUTO):
01775         plane_orientation = get_orientation(proxy)
01776 
01777     # Get cutting plane normal
01778     normal = None
01779 
01780     if (not is_planar_input(proxy)):
01781         normal = get_normal_by_orientation(plane_orientation,
01782                                            radians(angle1), radians(angle2))
01783 
01784         # Create slice filter
01785         slice_filter = pv.Slice(merge_blocks)
01786         slice_filter.SliceType = "Plane"
01787 
01788         # Set cutting plane normal
01789         slice_filter.SliceType.Normal = normal
01790 
01791         # Set cutting plane position
01792         if (is_relative):
01793             base_position = get_positions(1, normal,
01794                                           get_bounds(proxy), position)
01795             slice_filter.SliceOffsetValues = base_position
01796         else:
01797             slice_filter.SliceOffsetValues = position
01798 
01799         slice_filter.UpdatePipeline()
01800         poly_data = slice_filter
01801     else:
01802         normal = get_normal_by_orientation(plane_orientation, 0, 0)
01803 
01804     use_normal = 0
01805     # Geometry filter
01806     if not poly_data or poly_data.GetDataInformation().GetNumberOfCells() == 0:
01807         geometry_filter = pv.GeometryFilter(merge_blocks)
01808         poly_data = geometry_filter
01809         use_normal = 1  # TODO(MZN): workaround
01810 
01811     warp_scalar = None
01812     plot3d = None
01813     source = poly_data
01814 
01815     if is_data_on_cells(poly_data, field_name):
01816         # Cell data to point data
01817         cell_to_point = pv.CellDatatoPointData(poly_data)
01818         cell_to_point.PassCellData = 1
01819         source = cell_to_point
01820 
01821     scalars = ['POINTS', field_name]
01822 
01823     # Transform vector array to scalar array if necessary
01824     if (nb_components > 1):
01825         calc = get_calc_magnitude(source, EntityType.NODE, field_name)
01826         scalars = ['POINTS', calc.ResultArrayName]
01827         source = calc
01828 
01829     # Warp by scalar
01830     warp_scalar = pv.WarpByScalar(source)
01831     warp_scalar.Scalars = scalars
01832     warp_scalar.Normal = normal
01833     warp_scalar.UseNormal = use_normal
01834     if scale_factor is not None:
01835         warp_scalar.ScaleFactor = scale_factor
01836     else:
01837         def_scale = get_default_scale(PrsTypeEnum.PLOT3D,
01838                                       proxy, entity, field_name)
01839         warp_scalar.ScaleFactor = def_scale
01840 
01841     warp_scalar.UpdatePipeline()
01842     source = warp_scalar
01843 
01844     if (is_contour):
01845         # Contours
01846         contour = pv.Contour(warp_scalar)
01847         contour.PointMergeMethod = "Uniform Binning"
01848         contour.ContourBy = ['POINTS', field_name]
01849         scalar_range = get_data_range(proxy, entity,
01850                                       field_name, vector_mode)
01851         contour.Isosurfaces = get_contours(scalar_range, nb_contours)
01852         contour.UpdatePipeline()
01853         source = contour
01854 
01855     # Get Plot 3D representation object
01856     plot3d = pv.GetRepresentation(source)
01857 
01858     # Get lookup table
01859     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01860 
01861     # Set field range if necessary
01862     data_range = get_data_range(proxy, entity,
01863                                 field_name, vector_mode)
01864     lookup_table.LockScalarRange = 1
01865     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01866 
01867     # Set properties
01868     plot3d.ColorAttributeType = EntityType.get_pvtype(entity)
01869     plot3d.ColorArrayName = field_name
01870     plot3d.LookupTable = lookup_table
01871 
01872     # Add scalar bar
01873     add_scalar_bar(field_name, nb_components,
01874                    vector_mode, lookup_table, time_value)
01875 
01876     return plot3d
01877 
01878 
01879 def IsoSurfacesOnField(proxy, entity, field_name, timestamp_nb,
01880                        custom_range=None, nb_surfaces=10,
01881                        is_colored=True, color=None, vector_mode='Magnitude'):
01882     """Creates Iso Surfaces presentation on the given field.
01883 
01884     Arguments:
01885       proxy: the pipeline object, containig data
01886       entity: the entity type from PrsTypeEnum
01887       field_name: the field name
01888       timestamp_nb: the number of time step (1, 2, ...)
01889       custom_range: scalar range, if undefined the source range will be applied
01890       nb_surfaces: number of surfaces, which will be generated
01891       is_colored: this option allows to color the presentation according to
01892       the corresponding data array values. If False - the presentation will
01893       be one-coloured.
01894       color: defines the presentation color as [R, G, B] triple. Taken into
01895       account only if is_colored is False.
01896       vector_mode: the mode of transformation of vector values
01897       into scalar values, applicable only if the field contains vector values.
01898       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
01899 
01900     Returns:
01901       Iso Surfaces as representation object.
01902 
01903     """
01904     # We don't need mesh parts with no data on them
01905     if entity == EntityType.NODE:
01906         select_cells_with_data(proxy, on_points=[field_name])
01907     else:
01908         select_cells_with_data(proxy, on_cells=[field_name])
01909 
01910     # Check vector mode
01911     nb_components = get_nb_components(proxy, entity, field_name)
01912     check_vector_mode(vector_mode, nb_components)
01913 
01914     # Get time value
01915     time_value = get_time(proxy, timestamp_nb)
01916 
01917     # Set timestamp
01918     pv.GetRenderView().ViewTime = time_value
01919     pv.UpdatePipeline(time_value, proxy)
01920 
01921     # Extract only groups with data for the field
01922     new_proxy = extract_groups_for_field(proxy, field_name, entity)
01923 
01924     # Do merge
01925     source = pv.MergeBlocks(new_proxy)
01926 
01927     # Transform cell data into point data if necessary
01928     if is_data_on_cells(proxy, field_name):
01929         cell_to_point = pv.CellDatatoPointData(source)
01930         cell_to_point.PassCellData = 1
01931         source = cell_to_point
01932 
01933     contour_by = ['POINTS', field_name]
01934 
01935     # Transform vector array to scalar array if necessary
01936     if (nb_components > 1):
01937         calc = get_calc_magnitude(source, EntityType.NODE, field_name)
01938         contour_by = ['POINTS', calc.ResultArrayName]
01939         source = calc
01940 
01941     # Contour filter settings
01942     contour = pv.Contour(source)
01943     contour.ComputeScalars = 1
01944     contour.ContourBy = contour_by
01945 
01946     # Specify the range
01947     scalar_range = custom_range
01948     if (scalar_range is None):
01949         scalar_range = get_data_range(proxy, entity,
01950                                       field_name, cut_off=True)
01951 
01952     # Get contour values for the range
01953     surfaces = get_contours(scalar_range, nb_surfaces)
01954 
01955     # Set contour values
01956     contour.Isosurfaces = surfaces
01957 
01958     # Get Iso Surfaces representation object
01959     isosurfaces = pv.GetRepresentation(contour)
01960 
01961     # Get lookup table
01962     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
01963 
01964     # Set field range if necessary
01965     data_range = get_data_range(proxy, entity,
01966                                 field_name, vector_mode)
01967     lookup_table.LockScalarRange = 1
01968     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
01969 
01970     # Set display properties
01971     if (is_colored):
01972         isosurfaces.ColorAttributeType = EntityType.get_pvtype(entity)
01973         isosurfaces.ColorArrayName = field_name
01974     else:
01975         isosurfaces.ColorArrayName = ''
01976         if color:
01977             isosurfaces.DiffuseColor = color
01978     isosurfaces.LookupTable = lookup_table
01979 
01980     # Add scalar bar
01981     add_scalar_bar(field_name, nb_components,
01982                    vector_mode, lookup_table, time_value)
01983 
01984     return isosurfaces
01985 
01986 
01987 def GaussPointsOnField(proxy, entity, field_name,
01988                        timestamp_nb,
01989                        is_deformed=True, scale_factor=None,
01990                        is_colored=True, color=None,
01991                        primitive=GaussType.SPRITE,
01992                        is_proportional=True,
01993                        max_pixel_size=256,
01994                        multiplier=None, vector_mode='Magnitude'):
01995     """Creates Gauss Points on the given field.
01996 
01997     Arguments:
01998 
01999     proxy: the pipeline object, containig data
02000     entity: the field entity type from PrsTypeEnum
02001     field_name: the field name
02002     timestamp_nb: the number of time step (1, 2, ...)
02003     is_deformed: defines whether the Gauss Points will be deformed or not
02004     scale_factor -- the scale factor for deformation. Will be taken into
02005     account only if is_deformed is True.
02006     If not passed by user, default scale will be computed.
02007     is_colored -- defines whether the Gauss Points will be multicolored,
02008     using the corresponding data values
02009     color: defines the presentation color as [R, G, B] triple. Taken into
02010     account only if is_colored is False.
02011     primitive: primitive type from GaussType
02012     is_proportional: if True, the size of primitives will depends on
02013     the gauss point value
02014     max_pixel_size: the maximum sizr of the Gauss Points primitive in pixels
02015     multiplier: coefficient between data values and the size of primitives
02016     If not passed by user, default scale will be computed.
02017     vector_mode: the mode of transformation of vector values into
02018     scalar values, applicable only if the field contains vector values.
02019     Possible modes: 'Magnitude' - vector module;
02020     'X', 'Y', 'Z' - vector components.
02021 
02022     Returns:
02023       Gauss Points as representation object.
02024 
02025     """
02026     # We don't need mesh parts with no data on them
02027     if entity == EntityType.NODE:
02028         select_cells_with_data(proxy, on_points=[field_name])
02029     else:
02030         select_cells_with_data(proxy, on_cells=[field_name])
02031 
02032     # Check vector mode
02033     nb_components = get_nb_components(proxy, entity, field_name)
02034     check_vector_mode(vector_mode, nb_components)
02035 
02036     # Get time value
02037     time_value = get_time(proxy, timestamp_nb)
02038 
02039     # Set timestamp
02040     pv.GetRenderView().ViewTime = time_value
02041     proxy.UpdatePipeline(time=time_value)
02042 
02043     # Extract only groups with data for the field
02044     source = extract_groups_for_field(proxy, field_name, entity)
02045 
02046     # Quadrature point arrays
02047     qp_arrays = proxy.QuadraturePointArrays.Available
02048 
02049     # If no quadrature point array is passed, use cell centers
02050     if field_name in qp_arrays:
02051         generate_qp = pv.GenerateQuadraturePoints(source)
02052         generate_qp.SelectSourceArray = ['CELLS', 'ELGA_Offset']
02053         source = generate_qp
02054     else:
02055         # Cell centers
02056         cell_centers = pv.CellCenters(source)
02057         cell_centers.VertexCells = 1
02058         source = cell_centers
02059 
02060     source.UpdatePipeline()
02061 
02062     # Check if deformation enabled
02063     if is_deformed and nb_components > 1:
02064         vector_array = field_name
02065         # If the given vector array has only 2 components, add the third one
02066         if nb_components == 2:
02067             calc = get_add_component_calc(source,
02068                                           EntityType.NODE, field_name)
02069             vector_array = calc.ResultArrayName
02070             source = calc
02071 
02072         # Warp by vector
02073         warp_vector = pv.WarpByVector(source)
02074         warp_vector.Vectors = [vector_array]
02075         if scale_factor is not None:
02076             warp_vector.ScaleFactor = scale_factor
02077         else:
02078             def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE, proxy,
02079                                           entity, field_name)
02080             warp_vector.ScaleFactor = def_scale
02081         warp_vector.UpdatePipeline()
02082         source = warp_vector
02083 
02084     # Get Gauss Points representation object
02085     gausspnt = pv.GetRepresentation(source)
02086 
02087     # Get lookup table
02088     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
02089 
02090     # Set field range if necessary
02091     data_range = get_data_range(proxy, entity,
02092                                 field_name, vector_mode)
02093     lookup_table.LockScalarRange = 1
02094     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
02095 
02096     # Set display properties
02097     if is_colored:
02098         gausspnt.ColorAttributeType = EntityType.get_pvtype(entity)
02099         gausspnt.ColorArrayName = field_name
02100     else:
02101         gausspnt.ColorArrayName = ''
02102         if color:
02103             gausspnt.DiffuseColor = color
02104 
02105     gausspnt.LookupTable = lookup_table
02106 
02107     # Add scalar bar
02108     add_scalar_bar(field_name, nb_components,
02109                    vector_mode, lookup_table, time_value)
02110 
02111     # Set point sprite representation
02112     gausspnt.Representation = 'Point Sprite'
02113 
02114     # Point sprite settings
02115     gausspnt.InterpolateScalarsBeforeMapping = 0
02116     gausspnt.MaxPixelSize = max_pixel_size
02117 
02118     # Render mode
02119     gausspnt.RenderMode = GaussType.get_mode(primitive)
02120 
02121     #if primitive == GaussType.SPRITE:
02122         # Set texture
02123         # TODO(MZN): replace with pvsimple high-level interface
02124     #    texture = sm.CreateProxy("textures", "SpriteTexture")
02125     #    alphamprop = texture.GetProperty("AlphaMethod")
02126     #    alphamprop.SetElement(0, 2)  # Clamp
02127     #    alphatprop = texture.GetProperty("AlphaThreshold")
02128     #    alphatprop.SetElement(0, 63)
02129     #    maxprop = texture.GetProperty("Maximum")
02130     #    maxprop.SetElement(0, 255)
02131     #    texture.UpdateVTKObjects()
02132 
02133     #    gausspnt.Texture = texture
02134         #gausspnt.Texture.AlphaMethod = 'Clamp'
02135         #gausspnt.Texture.AlphaThreshold = 63
02136         #gausspnt.Texture.Maximum= 255
02137 
02138     # Proportional radius
02139     gausspnt.RadiusUseScalarRange = 0
02140     gausspnt.RadiusIsProportional = 0
02141 
02142     if is_proportional:
02143         mult = multiplier
02144         if mult is None:
02145             mult = abs(0.1 / data_range[1])
02146 
02147         gausspnt.RadiusScalarRange = data_range
02148         gausspnt.RadiusTransferFunctionEnabled = 1
02149         gausspnt.RadiusMode = 'Scalar'
02150         gausspnt.RadiusArray = ['POINTS', field_name]
02151         if nb_components > 1:
02152             v_comp = get_vector_component(vector_mode)
02153             gausspnt.RadiusVectorComponent = v_comp
02154         gausspnt.RadiusTransferFunctionMode = 'Table'
02155         gausspnt.RadiusScalarRange = data_range
02156         gausspnt.RadiusUseScalarRange = 1
02157         gausspnt.RadiusIsProportional = 1
02158         gausspnt.RadiusProportionalFactor = mult
02159     else:
02160         gausspnt.RadiusTransferFunctionEnabled = 0
02161         gausspnt.RadiusMode = 'Constant'
02162         gausspnt.RadiusArray = ['POINTS', 'Constant Radius']
02163 
02164     return gausspnt
02165 
02166 
02167 def StreamLinesOnField(proxy, entity, field_name, timestamp_nb,
02168                        direction='BOTH', is_colored=False, color=None,
02169                        vector_mode='Magnitude'):
02170     """Creates Stream Lines presentation on the given field.
02171 
02172     Arguments:
02173       proxy: the pipeline object, containig data
02174       entity: the entity type from PrsTypeEnum
02175       field_name: the field name
02176       timestamp_nb: the number of time step (1, 2, ...)
02177       direction: the stream lines direction ('FORWARD', 'BACKWARD' or 'BOTH')
02178       is_colored: this option allows to color the presentation according to
02179       the corresponding data values. If False - the presentation will
02180       be one-coloured.
02181       color: defines the presentation color as [R, G, B] triple. Taken into
02182       account only if is_colored is False.
02183       vector_mode: the mode of transformation of vector values
02184       into scalar values, applicable only if the field contains vector values.
02185       Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
02186 
02187     Returns:
02188       Stream Lines as representation object.
02189 
02190     """
02191     # We don't need mesh parts with no data on them
02192     if entity == EntityType.NODE:
02193         select_cells_with_data(proxy, on_points=[field_name])
02194     else:
02195         select_cells_with_data(proxy, on_cells=[field_name])
02196 
02197     # Check vector mode
02198     nb_components = get_nb_components(proxy, entity, field_name)
02199     check_vector_mode(vector_mode, nb_components)
02200 
02201     # Get time value
02202     time_value = get_time(proxy, timestamp_nb)
02203 
02204     # Set timestamp
02205     pv.GetRenderView().ViewTime = time_value
02206     pv.UpdatePipeline(time_value, proxy)
02207 
02208     # Extract only groups with data for the field
02209     new_proxy = extract_groups_for_field(proxy, field_name, entity)
02210 
02211     # Do merge
02212     source = pv.MergeBlocks(new_proxy)
02213 
02214     # Cell data to point data
02215     if is_data_on_cells(proxy, field_name):
02216         cell_to_point = pv.CellDatatoPointData(source)
02217         cell_to_point.PassCellData = 1
02218         cell_to_point.UpdatePipeline()
02219         source = cell_to_point
02220 
02221     vector_array = field_name
02222     # If the given vector array has only 2 components, add the third one
02223     if nb_components == 2:
02224         calc = get_add_component_calc(source, EntityType.NODE, field_name)
02225         vector_array = calc.ResultArrayName
02226         calc.UpdatePipeline()
02227         source = calc
02228 
02229     # Stream Tracer
02230     stream = pv.StreamTracer(source)
02231     stream.SeedType = "Point Source"
02232     stream.Vectors = ['POINTS', vector_array]
02233     stream.SeedType = "Point Source"
02234     stream.IntegrationDirection = direction
02235     stream.IntegratorType = 'Runge-Kutta 2'
02236     stream.UpdatePipeline()
02237 
02238     # Get Stream Lines representation object
02239     if is_empty(stream):
02240         return None
02241     streamlines = pv.GetRepresentation(stream)
02242 
02243     # Get lookup table
02244     lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
02245 
02246     # Set field range if necessary
02247     data_range = get_data_range(new_proxy, entity,
02248                                 field_name, vector_mode)
02249     lookup_table.LockScalarRange = 1
02250     lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
02251 
02252     # Set properties
02253     if is_colored:
02254         streamlines.ColorAttributeType = EntityType.get_pvtype(entity)
02255         streamlines.ColorArrayName = field_name
02256     else:
02257         streamlines.ColorArrayName = ''
02258         if color:
02259             streamlines.DiffuseColor = color
02260 
02261     streamlines.LookupTable = lookup_table
02262 
02263     # Add scalar bar
02264     add_scalar_bar(field_name, nb_components,
02265                    vector_mode, lookup_table, time_value)
02266 
02267     return streamlines
02268 
02269 
02270 def MeshOnEntity(proxy, mesh_name, entity):
02271     """Creates submesh of the entity type for the mesh.
02272 
02273     Arguments:
02274       proxy -- the pipeline object, containig data
02275       mesh_name -- the mesh name
02276       entity -- the entity type
02277 
02278     Returns:
02279       Submesh as representation object of the given source.
02280 
02281     """
02282     # Select all cell types
02283     select_all_cells(proxy)
02284 
02285     # Get subset of groups on the given entity
02286     subset = get_group_names(proxy, mesh_name, entity)
02287 
02288     # Select only groups of the given entity type
02289     proxy.Groups = subset
02290     proxy.UpdatePipeline()
02291 
02292     # Get representation object if the submesh is not empty
02293     prs = None
02294     if (proxy.GetDataInformation().GetNumberOfPoints() or
02295         proxy.GetDataInformation().GetNumberOfCells()):
02296         prs = pv.GetRepresentation(proxy)
02297         prs.ColorArrayName = ''
02298 
02299     return prs
02300 
02301 
02302 def MeshOnGroup(proxy, group_name):
02303     """Creates submesh on the group.
02304 
02305     Arguments:
02306       proxy -- the pipeline object, containig data
02307       group_name -- the full group name
02308 
02309     Returns:
02310       Representation object of the given source with single group
02311       selected.
02312 
02313     """
02314     # Select all cell types
02315     select_all_cells(proxy)
02316 
02317     # Select only the group with the given name
02318     one_group = [group_name]
02319     proxy.Groups = one_group
02320     proxy.UpdatePipeline()
02321 
02322     # Get representation object if the submesh is not empty
02323     prs = None
02324 
02325     # Check if the group was set
02326     if proxy.Groups.GetData() == one_group:
02327         group_entity = get_group_entity(group_name)
02328         # Check if the submesh is not empty
02329         nb_items = 0
02330         if group_entity == EntityType.NODE:
02331             nb_items = proxy.GetDataInformation().GetNumberOfPoints()
02332         elif group_entity == EntityType.CELL:
02333             nb_items = proxy.GetDataInformation().GetNumberOfCells()
02334 
02335         if nb_items:
02336             prs = pv.GetRepresentation(proxy)
02337             prs.ColorArrayName = ''
02338 
02339     return prs
02340 
02341 
02342 def CreatePrsForFile(paravis_instance, file_name, prs_types,
02343                      picture_dir, picture_ext):
02344     """Build presentations of the given types for the file.
02345 
02346     Build presentations for all fields on all timestamps.
02347 
02348     Arguments:
02349       paravis_instance: ParaVis module instance object
02350       file_name: full path to the MED file
02351       prs_types: the list of presentation types to build
02352       picture_dir: the directory path for saving snapshots
02353       picture_ext: graphics files extension (determines file type)
02354 
02355     """
02356     # Import MED file
02357     print "Import " + file_name.split(os.sep)[-1] + "..."
02358 
02359     try:
02360         paravis_instance.ImportFile(file_name)
02361         proxy = pv.GetActiveSource()
02362         if proxy is None:
02363             print "FAILED"
02364         else:
02365             proxy.UpdatePipeline()
02366             print "OK"
02367     except:
02368         print "FAILED"
02369     else:
02370         # Get view
02371         view = pv.GetRenderView()
02372 
02373         # Create required presentations for the proxy
02374         CreatePrsForProxy(proxy, view, prs_types,
02375                           picture_dir, picture_ext)
02376 
02377 
02378 def CreatePrsForProxy(proxy, view, prs_types, picture_dir, picture_ext):
02379     """Build presentations of the given types for all fields of the proxy.
02380 
02381     Save snapshots in graphics files (type depends on the given extension).
02382     Stores the files in the given directory.
02383 
02384     Arguments:
02385       proxy: the pipeline object, containig data
02386       view: the render view
02387       prs_types: the list of presentation types to build
02388       picture_dir: the directory path for saving snapshots
02389       picture_ext: graphics files extension (determines file type)
02390 
02391     """
02392     # List of the field names
02393     field_names = list(proxy.PointArrays.GetData())
02394     nb_on_nodes = len(field_names)
02395     field_names.extend(proxy.CellArrays.GetData())
02396 
02397     # Add path separator to the end of picture path if necessery
02398     if not picture_dir.endswith(os.sep):
02399         picture_dir += os.sep
02400 
02401     # Mesh Presentation
02402     if PrsTypeEnum.MESH in prs_types:
02403         # Create Mesh presentation. Build all possible submeshes.
02404 
02405         # Remember the current state
02406         groups = list(proxy.Groups)
02407 
02408         # Iterate on meshes
02409         mesh_names = get_mesh_names(proxy)
02410         for mesh_name in mesh_names:
02411             # Build mesh on nodes and cells
02412             for entity in (EntityType.NODE, EntityType.CELL):
02413                 entity_name = EntityType.get_name(entity)
02414                 if if_possible(proxy, mesh_name, entity, PrsTypeEnum.MESH):
02415                     print "Creating submesh on " + entity_name + " for '" + mesh_name + "' mesh... "
02416                     prs = MeshOnEntity(proxy, mesh_name, entity)
02417                     if prs is None:
02418                         print "FAILED"
02419                         continue
02420                     else:
02421                         print "OK"
02422                     # Construct image file name
02423                     pic_name = picture_dir + mesh_name + "_" + entity_name + "." + picture_ext
02424 
02425                     # Show and dump the presentation into a graphics file
02426                     process_prs_for_test(prs, view, pic_name, False)
02427 
02428                 # Build submesh on all groups of the mesh
02429                 mesh_groups = get_group_names(proxy, mesh_name,
02430                                               entity, wo_nogroups=True)
02431                 for group in mesh_groups:
02432                     print "Creating submesh on group " + group + "... "
02433                     prs = MeshOnGroup(proxy, group)
02434                     if prs is None:
02435                         print "FAILED"
02436                         continue
02437                     else:
02438                         print "OK"
02439                     # Construct image file name
02440                     pic_name = picture_dir + group.replace('/', '_') + "." + picture_ext
02441 
02442                     # Show and dump the presentation into a graphics file
02443                     process_prs_for_test(prs, view, pic_name, False)
02444 
02445         # Restore the state
02446         proxy.Groups = groups
02447         proxy.UpdatePipeline()
02448 
02449     # Presentations on fields
02450     for (i, field_name) in enumerate(field_names):
02451         # Select only the current field:
02452         # necessary for getting the right timestamps
02453         cell_arrays = proxy.CellArrays.GetData()
02454         point_arrays = proxy.PointArrays.GetData()
02455         field_entity = None
02456         if (i >= nb_on_nodes):
02457             field_entity = EntityType.CELL
02458             proxy.PointArrays.DeselectAll()
02459             proxy.CellArrays = [field_name]
02460         else:
02461             field_entity = EntityType.NODE
02462             proxy.CellArrays.DeselectAll()
02463             proxy.PointArrays = [field_name]
02464 
02465         # Get timestamps
02466         proxy.UpdatePipelineInformation()
02467         timestamps = proxy.TimestepValues.GetData()
02468 
02469         # Restore fields selection state
02470         proxy.CellArrays = cell_arrays
02471         proxy.PointArrays = point_arrays
02472         proxy.UpdatePipelineInformation()
02473 
02474         for prs_type in prs_types:
02475             # Ignore mesh presentation
02476             if prs_type == PrsTypeEnum.MESH:
02477                 continue
02478 
02479             # Get name of presentation type
02480             prs_name = PrsTypeEnum.get_name(prs_type)
02481 
02482             # Build the presentation if possible
02483             possible = if_possible(proxy, field_name,
02484                                    field_entity, prs_type)
02485             if possible:
02486                 # Presentation type for graphics file name
02487                 f_prs_type = prs_name.replace(' ', '').upper()
02488 
02489                 for timestamp_nb in xrange(1, len(timestamps) + 1):
02490                     time = timestamps[timestamp_nb - 1]
02491                     print "Creating " + prs_name + " on " + field_name + ", time = " + str(time) + "... "
02492                     prs = create_prs(prs_type, proxy,
02493                                      field_entity, field_name, timestamp_nb)
02494                     if prs is None:
02495                         print "FAILED"
02496                         continue
02497                     else:
02498                         print "OK"
02499 
02500                     # Construct image file name
02501                     pic_name = picture_dir + field_name + "_" + str(time) + "_" + f_prs_type + "." + picture_ext
02502 
02503                     # Show and dump the presentation into a graphics file
02504                     process_prs_for_test(prs, view, pic_name)