1 # Copyright (C) 2010-2015 CEA/DEN, EDF R&D
3 # This library is free software; you can redistribute it and/or
4 # modify it under the terms of the GNU Lesser General Public
5 # License as published by the Free Software Foundation; either
6 # version 2.1 of the License, or (at your option) any later version.
8 # This library is distributed in the hope that it will be useful,
9 # but WITHOUT ANY WARRANTY; without even the implied warranty of
10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 # Lesser General Public License for more details.
13 # You should have received a copy of the GNU Lesser General Public
14 # License along with this library; if not, write to the Free Software
15 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 This module is intended to provide Python API for building presentations
22 typical for Post-Pro module (Scalar Map, Deformed Shape, Vectors, etc.)
26 from __future__ import division
27 ##from __future__ import print_function
32 from math import sqrt, sin, cos, radians
33 from string import upper
35 # Do not use pv as a short name.
36 # It is a name of function from numpy and may be redefined implicitly by 'from numpy import *' call.
37 # import pvsimple as pv
38 import pvsimple as pvs
40 # # TODO(MZN): to be removed (issue with Point Sprite texture)
41 # #import paravisSM as sm
43 # import paraview.simple as pvs
44 # import paraview.servermanager as sm
50 VTK_LARGE_FLOAT = 1E+38
51 GAP_COEFFICIENT = 0.0001
56 _med_field_sep = '@@][@@'
62 Post-Pro presentation types.
70 DEFORMEDSHAPESCALARMAP = 6
76 _type2name = {MESH: 'Mesh',
77 SCALARMAP: 'Scalar Map',
78 ISOSURFACES: 'Iso Surfaces',
79 CUTPLANES: 'Cut Planes',
80 CUTLINES: 'Cut Lines',
81 DEFORMEDSHAPE: 'Deformed Shape',
82 DEFORMEDSHAPESCALARMAP: 'Deformed Shape And Scalar Map',
85 STREAMLINES: 'Stream Lines',
86 GAUSSPOINTS: 'Gauss Points'}
89 def get_name(cls, type):
90 """Return presentaion name by its type."""
91 return cls._type2name[type]
101 _type2name = {NODE: 'P1',
104 _name2type = {'P1': NODE,
107 _type2pvtype = {NODE: 'POINT_DATA',
111 def get_name(cls, type):
112 """Return entity name (used in full group names) by its type."""
113 return cls._type2name[type]
116 def get_type(cls, name):
117 """Return entity type by its name (used in full group names)."""
118 return cls._name2type[name]
121 def get_pvtype(cls, type):
122 """Return entity type from ['CELL_DATA', 'POINT_DATA']"""
123 return cls._type2pvtype[type]
127 """Orientation types.
129 Defines a set of plane orientation possibilities:
130 AUTO: plane orientation should be calculated.
131 XY: plane formed by X and Y axis.
132 YZ: plane formed by Y and Z axis.
133 ZX: plane formed by Z and X axis
145 Set of elements defining the position of the vector head:
146 CENTER: in the center of the vector
147 TAIL: in the tail of the vector
148 HEAD: in the head of the vector
158 Gauss Points primitive types.
164 _type2mode = {SPRITE: 'Texture',
165 POINT: 'SimplePoint',
166 SPHERE: 'Sphere (Texture)'}
169 def get_mode(cls, type):
170 """Return paraview point sprite mode by the primitive type."""
171 return cls._type2mode[type]
174 # Auxiliary functions
176 def get_field_mesh_name(full_field_name):
177 """Return mesh name of the field by its full name."""
178 aList = full_field_name.split('/')
180 field_name = full_field_name.split('/')[1]
184 def get_field_entity(full_field_name):
185 """Return entity type of the field by its full name."""
186 aList = full_field_name.split(_med_field_sep)
188 entity_name = full_field_name.split(_med_field_sep)[-1]
189 entity = EntityType.get_type(entity_name)
193 def get_field_short_name(full_field_name):
194 """Return short name of the field by its full name."""
195 aList = full_field_name.split('/')
197 short_name_with_type = full_field_name.split('/')[-1]
198 short_name = short_name_with_type.split(_med_field_sep)[0]
202 def find_mesh_full_name(proxy, short_mesh_name):
203 """Return full mesh path by short mesh name, if found"""
204 proxy.UpdatePipeline()
205 all_mesh_names = get_mesh_full_names(proxy)
206 for name in all_mesh_names:
207 if short_mesh_name == get_field_short_name(name):
211 def process_prs_for_test(prs, view, picture_name, show_bar=True):
212 """Show presentation and record snapshot image.
215 prs: the presentation to show
216 view: the render view
217 picture_name: the full name of the graphics file to save
218 show_bar: to show scalar bar or not
221 # Show the presentation only
222 display_only(prs, view)
226 if show_bar and _current_bar:
227 _current_bar.Visibility = 1
232 # Create a directory for screenshot if necessary
233 file_name = re.sub("\s+", "_", picture_name)
234 pic_dir = os.path.dirname(picture_name)
235 if not os.path.exists(pic_dir):
239 print "Write image:", file_name
240 pvs.WriteImage(file_name, view=view, Magnification=1)
243 def reset_view(view=None):
246 Set predefined (taken from Post-Pro) camera settings.
247 If the view is not passed, the active view is used.
251 view = pvs.GetRenderView()
254 view.CameraFocalPoint = [0.0, 0.0, 0.0]
255 view.CameraViewUp = [0.0, 0.0, 1.0]
256 view.CameraPosition = [738.946, -738.946, 738.946]
258 # Turn on the headligth
260 view.LightIntensity = 0.5
262 # Use parallel projection
263 view.CameraParallelProjection = 1
266 pvs.Render(view=view)
269 def hide_all(view, to_remove=False):
270 """Hide all representations in the view."""
272 view = pvs.GetRenderView()
274 rep_list = view.Representations
276 if hasattr(rep, 'Visibility') and rep.Visibility != 0:
279 view.Representations.remove(rep)
280 pvs.Render(view=view)
283 def display_only(prs, view=None):
284 """Display only the given presentation in the view."""
286 view = pvs.GetRenderView()
288 rep_list = view.Representations
290 if hasattr(rep, 'Visibility'):
291 rep.Visibility = (rep == prs)
292 pvs.Render(view=view)
295 def set_visible_lines(xy_prs, lines):
296 """Set visible only the given lines for XYChartRepresentation."""
297 sv = xy_prs.GetProperty("SeriesVisibility").GetData()
300 for i in xrange(0, len(sv)):
303 if line_name in lines:
310 xy_prs.SeriesVisibility = sv
313 def check_vector_mode(vector_mode, nb_components):
314 """Check vector mode.
316 Check if vector mode is correct for the data array with the
317 given number of components.
320 vector_mode: 'Magnitude', 'X', 'Y' or 'Z'
321 nb_components: number of component in the data array
324 ValueError: in case of the vector mode is unexistent
328 if vector_mode not in ('Magnitude', 'X', 'Y', 'Z'):
329 raise ValueError("Unexistent vector mode: " + vector_mode)
331 if ((nb_components == 1 and (vector_mode == 'Y' or vector_mode == 'Z')) or
332 (nb_components == 2 and vector_mode == 'Z')):
333 raise ValueError("Incorrect vector mode " + vector_mode + " for " +
334 nb_components + "-component field")
337 def get_vector_component(vector_mode):
338 """Get vector component as ineger.
340 Translate vector component notation from string
350 if vector_mode == 'X':
352 elif vector_mode == 'Y':
354 elif vector_mode == 'Z':
360 def get_data_range(proxy, entity, field_name, vector_mode='Magnitude',
362 """Get data range for the field.
365 proxy: the pipeline object, containig data array for the field
366 entity: the field entity
367 field_name: the field name
368 vector_mode: the vector mode ('Magnitude', 'X', 'Y' or 'Z')
371 Data range as [min, max]
374 proxy.UpdatePipeline()
375 entity_data_info = None
376 field_data = proxy.GetFieldDataInformation()
378 if field_name in field_data.keys():
379 entity_data_info = field_data
380 elif entity == EntityType.CELL:
381 entity_data_info = proxy.GetCellDataInformation()
382 elif entity == EntityType.NODE:
383 entity_data_info = proxy.GetPointDataInformation()
387 if field_name in entity_data_info.keys():
388 vcomp = get_vector_component(vector_mode)
389 data_range = entity_data_info[field_name].GetComponentRange(vcomp)
391 pv_entity = EntityType.get_pvtype(entity)
392 warnings.warn("Field " + field_name +
393 " is unknown for " + pv_entity + "!")
396 if cut_off and (data_range[0] <= data_range[1]):
397 data_range = list(data_range)
398 delta = abs(data_range[1] - data_range[0]) * GAP_COEFFICIENT
399 data_range[0] += delta
400 data_range[1] -= delta
405 def get_bounds(proxy):
406 """Get bounds of the proxy in 3D."""
407 proxy.UpdatePipeline()
408 dataInfo = proxy.GetDataInformation()
409 bounds_info = dataInfo.GetBounds()
413 def get_x_range(proxy):
414 """Get X range of the proxy bounds in 3D."""
415 proxy.UpdatePipeline()
416 bounds_info = get_bounds(proxy)
417 return bounds_info[0:2]
420 def get_y_range(proxy):
421 """Get Y range of the proxy bounds in 3D."""
422 proxy.UpdatePipeline()
423 bounds_info = get_bounds(proxy)
424 return bounds_info[2:4]
427 def get_z_range(proxy):
428 """Get Z range of the proxy bounds in 3D."""
429 proxy.UpdatePipeline()
430 bounds_info = get_bounds(proxy)
431 return bounds_info[4:6]
434 def is_planar_input(proxy):
435 """Check if the given input is planar."""
436 proxy.UpdatePipeline()
437 bounds_info = get_bounds(proxy)
439 if (abs(bounds_info[0] - bounds_info[1]) <= FLT_MIN or
440 abs(bounds_info[2] - bounds_info[3]) <= FLT_MIN or
441 abs(bounds_info[4] - bounds_info[5]) <= FLT_MIN):
447 def is_data_on_cells(proxy, field_name):
448 """Check the existence of a field on cells with the given name."""
449 proxy.UpdatePipeline()
450 cell_data_info = proxy.GetCellDataInformation()
451 return (field_name in cell_data_info.keys())
455 """Check if the object contains any points or cells.
458 True: if the given proxy doesn't contain any points or cells
462 proxy.UpdatePipeline()
463 data_info = proxy.GetDataInformation()
465 nb_cells = data_info.GetNumberOfCells()
466 nb_points = data_info.GetNumberOfPoints()
468 return not(nb_cells + nb_points)
471 def get_orientation(proxy):
472 """Get the optimum cutting plane orientation for Plot 3D."""
473 proxy.UpdatePipeline()
474 orientation = Orientation.XY
476 bounds = get_bounds(proxy)
477 delta = [bounds[1] - bounds[0],
478 bounds[3] - bounds[2],
479 bounds[5] - bounds[4]]
481 if (delta[0] >= delta[1] and delta[0] >= delta[2]):
482 if (delta[1] >= delta[2]):
483 orientation = Orientation.XY
485 orientation = Orientation.ZX
486 elif (delta[1] >= delta[0] and delta[1] >= delta[2]):
487 if (delta[0] >= delta[2]):
488 orientation = Orientation.XY
490 orientation = Orientation.YZ
491 elif (delta[2] >= delta[0] and delta[2] >= delta[1]):
492 if (delta[0] >= delta[1]):
493 orientation = Orientation.ZX
495 orientation = Orientation.YZ
500 def dot_product(a, b):
501 """Dot product of two 3-vectors."""
502 dot = a[0] * b[0] + a[1] * b[1] + a[2] * b[2]
506 def multiply3x3(a, b):
507 """Mutltiply one 3x3 matrix by another."""
513 c[0][i] = a[0][0] * b[0][i] + a[0][1] * b[1][i] + a[0][2] * b[2][i]
514 c[1][i] = a[1][0] * b[0][i] + a[1][1] * b[1][i] + a[1][2] * b[2][i]
515 c[2][i] = a[2][0] * b[0][i] + a[2][1] * b[1][i] + a[2][2] * b[2][i]
521 """Get X rotation matrix by angle."""
522 rx = [[1.0, 0.0, 0.0],
523 [0.0, cos(ang), -sin(ang)],
524 [0.0, sin(ang), cos(ang)]]
530 """Get Y rotation matrix by angle."""
531 ry = [[cos(ang), 0.0, sin(ang)],
533 [-sin(ang), 0.0, cos(ang)]]
539 """Get Z rotation matrix by angle."""
540 rz = [[cos(ang), -sin(ang), 0.0],
541 [sin(ang), cos(ang), 0.0],
547 def get_normal_by_orientation(orientation, ang1=0, ang2=0):
548 """Get normal for the plane by its orientation."""
550 rotation = [[], [], []]
551 rx = ry = rz = [[1.0, 0.0, 0.0],
555 normal = [0.0, 0.0, 0.0]
556 if orientation == Orientation.XY:
561 rotation = multiply3x3(rx, ry)
563 elif orientation == Orientation.ZX:
568 rotation = multiply3x3(rz, rx)
570 elif orientation == Orientation.YZ:
575 rotation = multiply3x3(ry, rz)
578 for i in xrange(0, 3):
579 normal[i] = rotation[i][i_plane]
584 def get_bound_project(bound_box, dir):
585 """Get bounds projection"""
586 bound_points = [[bound_box[0], bound_box[2], bound_box[4]],
587 [bound_box[1], bound_box[2], bound_box[4]],
588 [bound_box[0], bound_box[3], bound_box[4]],
589 [bound_box[1], bound_box[3], bound_box[4]],
590 [bound_box[0], bound_box[2], bound_box[5]],
591 [bound_box[1], bound_box[2], bound_box[5]],
592 [bound_box[0], bound_box[3], bound_box[5]],
593 [bound_box[1], bound_box[3], bound_box[5]]]
595 bound_prj = [0, 0, 0]
596 bound_prj[0] = dot_product(dir, bound_points[0])
597 bound_prj[1] = bound_prj[0]
599 for i in xrange(1, 8):
600 tmp = dot_product(dir, bound_points[i])
601 if bound_prj[1] < tmp:
603 if bound_prj[0] > tmp:
606 bound_prj[2] = bound_prj[1] - bound_prj[0]
607 bound_prj[1] = bound_prj[0] + (1.0 - EPS) * bound_prj[2]
608 bound_prj[0] = bound_prj[0] + EPS * bound_prj[2]
609 bound_prj[2] = bound_prj[1] - bound_prj[0]
614 def get_positions(nb_planes, dir, bounds, displacement):
615 """Compute plane positions."""
617 bound_prj = get_bound_project(bounds, dir)
619 step = bound_prj[2] / (nb_planes - 1)
620 abs_displacement = step * displacement
621 start_pos = bound_prj[0] - 0.5 * step + abs_displacement
622 for i in xrange(nb_planes):
623 pos = start_pos + i * step
624 positions.append(pos)
626 pos = bound_prj[0] + bound_prj[2] * displacement
627 positions.append(pos)
632 def get_contours(scalar_range, nb_contours):
633 """Generate contour values."""
635 for i in xrange(nb_contours):
636 pos = scalar_range[0] + i * (
637 scalar_range[1] - scalar_range[0]) / (nb_contours - 1)
643 def get_nb_components(proxy, entity, field_name):
644 """Return number of components for the field."""
645 proxy.UpdatePipeline()
646 entity_data_info = None
647 field_data = proxy.GetFieldDataInformation()
649 if field_name in field_data.keys():
650 entity_data_info = field_data
651 elif entity == EntityType.CELL:
652 select_cells_with_data(proxy, on_cells=[field_name])
653 entity_data_info = proxy.GetCellDataInformation()
654 elif entity == EntityType.NODE:
655 select_cells_with_data(proxy, on_points=[field_name])
656 entity_data_info = proxy.GetPointDataInformation()
659 if field_name in entity_data_info.keys():
660 nb_comp = entity_data_info[field_name].GetNumberOfComponents()
662 pv_entity = EntityType.get_pvtype(entity)
663 raise ValueError("Field " + field_name +
664 " is unknown for " + pv_entity + "!")
669 def get_scale_factor(proxy):
670 """Compute scale factor."""
674 proxy.UpdatePipeline()
675 data_info = proxy.GetDataInformation()
677 nb_cells = data_info.GetNumberOfCells()
678 nb_points = data_info.GetNumberOfPoints()
679 nb_elements = nb_cells if nb_cells > 0 else nb_points
680 bounds = get_bounds(proxy)
685 for i in xrange(0, 6, 2):
686 vol = abs(bounds[i + 1] - bounds[i])
691 if nb_elements == 0 or dim < 1 / VTK_LARGE_FLOAT:
694 volume /= nb_elements
696 return pow(volume, 1 / dim)
699 def get_default_scale(prs_type, proxy, entity, field_name):
700 """Get default scale factor."""
701 proxy.UpdatePipeline()
702 data_range = get_data_range(proxy, entity, field_name)
704 if prs_type == PrsTypeEnum.DEFORMEDSHAPE:
705 EPS = 1.0 / VTK_LARGE_FLOAT
706 if abs(data_range[1]) > EPS:
707 scale_factor = get_scale_factor(proxy)
708 return scale_factor / data_range[1]
709 elif prs_type == PrsTypeEnum.PLOT3D:
710 bounds = get_bounds(proxy)
711 length = sqrt((bounds[1] - bounds[0]) ** 2 +
712 (bounds[3] - bounds[2]) ** 2 +
713 (bounds[5] - bounds[4]) ** 2)
716 if data_range[1] > 0:
717 return length / data_range[1] * EPS
722 def get_calc_magnitude(proxy, array_entity, array_name):
723 """Compute magnitude for the given vector array via Calculator.
726 the calculator object.
729 proxy.UpdatePipeline()
732 # Transform vector array to scalar array if possible
733 nb_components = get_nb_components(proxy, array_entity, array_name)
734 if (nb_components > 1):
735 calculator = pvs.Calculator(proxy)
736 attribute_mode = "Point Data"
737 if array_entity != EntityType.NODE:
738 attribute_mode = "Cell Data"
739 calculator.AttributeMode = attribute_mode
740 if (nb_components == 2):
741 # Workaroud: calculator unable to compute magnitude
742 # if number of components equal to 2
743 func = "sqrt(" + array_name + "_X^2+" + array_name + "_Y^2)"
744 calculator.Function = func
746 calculator.Function = "mag(" + array_name + ")"
747 calculator.ResultArrayName = array_name + "_magnitude"
748 calculator.UpdatePipeline()
753 def get_add_component_calc(proxy, array_entity, array_name):
754 """Creates 3-component array from 2-component.
756 The first two components is from the original array. The 3rd component
758 If the number of components is not equal to 2 - return original array name.
761 the calculator object.
764 proxy.UpdatePipeline()
767 nb_components = get_nb_components(proxy, array_entity, array_name)
768 if nb_components == 2:
769 calculator = pvs.Calculator(proxy)
770 attribute_mode = "Point Data"
771 if array_entity != EntityType.NODE:
772 attribute_mode = "Cell Data"
773 calculator.AttributeMode = attribute_mode
774 expression = "iHat * " + array_name + "_X + jHat * " + array_name + "_Y + kHat * 0"
775 calculator.Function = expression
776 calculator.ResultArrayName = array_name + "_3c"
777 calculator.UpdatePipeline()
782 def select_all_cells(proxy):
783 """Select all cell types.
785 Used in creation of mesh/submesh presentation.
788 proxy.UpdatePipeline()
789 extractCT = pvs.ExtractCellType()
790 extractCT.AllGeoTypes = extractCT.GetProperty("GeoTypesInfo")[::2]
791 extractCT.UpdatePipelineInformation()
794 def select_cells_with_data(proxy, on_points=[], on_cells=[], on_gauss=[]):
795 """Select cell types with data.
797 Only cell types with data for the given fields will be selected.
798 If no fields defined (neither on points nor on cells) only cell
799 types with data for even one field (from available) will be selected.
802 if not proxy.GetProperty("FieldsTreeInfo"):
805 proxy.UpdatePipeline()
806 if not hasattr(proxy, 'Entity'):
807 fields_info = proxy.GetProperty("FieldsTreeInfo")[::2]
808 arr_name_with_dis=[elt.split("/")[-1] for elt in fields_info]
811 for name in on_gauss:
812 fields.append(name+_med_field_sep+'GAUSS')
813 for name in on_cells:
814 fields.append(name+_med_field_sep+'P0')
815 for name in on_points:
816 fields.append(name+_med_field_sep+'P1')
820 if arr_name_with_dis.count(name) > 0:
821 index = arr_name_with_dis.index(name)
822 field_list.append(fields_info[index])
825 proxy.AllArrays = field_list
826 proxy.UpdatePipeline()
827 return len(field_list) != 0
829 # TODO: VTN. Looks like this code is out of date.
831 #all_cell_types = proxy.CellTypes.Available
832 all_cell_types = proxy.Entity.Available
833 all_arrays = list(proxy.CellArrays.GetData())
834 all_arrays.extend(proxy.PointArrays.GetData())
837 file_name = proxy.FileName.split(os.sep)[-1]
838 print "Warning: " + file_name + " doesn't contain any data array."
840 # List of cell types to be selected
843 for cell_type in all_cell_types:
844 #proxy.CellTypes = [cell_type]
845 proxy.Entity = [cell_type]
846 proxy.UpdatePipeline()
848 cell_arrays = proxy.GetCellDataInformation().keys()
849 point_arrays = proxy.GetPointDataInformation().keys()
851 if on_points or on_cells:
852 if on_points is None:
857 if (all(array in cell_arrays for array in on_cells) and
858 all(array in point_arrays for array in on_points)):
859 # Add cell type to the list
860 cell_types_on.append(cell_type)
862 in_arrays = lambda array: ((array in cell_arrays) or
863 (array in point_arrays))
864 if any(in_arrays(array) for array in all_arrays):
865 cell_types_on.append(cell_type)
868 #proxy.CellTypes = cell_types_on
869 proxy.Entity = cell_types_on
870 proxy.UpdatePipeline()
872 def if_possible(proxy, field_name, entity, prs_type, extrGrps=None):
873 """Check if the presentation creation is possible on the given field."""
874 proxy.UpdatePipeline()
876 if (prs_type == PrsTypeEnum.DEFORMEDSHAPE or
877 prs_type == PrsTypeEnum.DEFORMEDSHAPESCALARMAP or
878 prs_type == PrsTypeEnum.VECTORS or
879 prs_type == PrsTypeEnum.STREAMLINES):
880 nb_comp = get_nb_components(proxy, entity, field_name)
881 result = (nb_comp > 1)
882 elif (prs_type == PrsTypeEnum.GAUSSPOINTS):
883 result = (entity == EntityType.CELL or
884 field_name in proxy.QuadraturePointArrays.Available)
885 elif (prs_type == PrsTypeEnum.MESH):
886 result = len(get_group_names(extrGrps)) > 0
891 def add_scalar_bar(field_name, nb_components,
892 vector_mode, lookup_table, time_value):
893 """Add scalar bar with predefined properties."""
896 # Construct bar title
897 title = "\n".join([field_name, str(time_value)])
898 if nb_components > 1:
899 title = "\n".join([title, vector_mode])
902 scalar_bar = pvs.CreateScalarBar(Enabled=1)
903 scalar_bar.Orientation = 'Vertical'
904 scalar_bar.Title = title
905 scalar_bar.LookupTable = lookup_table
907 # Set default properties same as in Post-Pro
908 scalar_bar.NumberOfLabels = 5
909 scalar_bar.AutomaticLabelFormat = 0
910 scalar_bar.LabelFormat = '%-#6.6g'
912 scalar_bar.TitleFontFamily = 'Arial'
913 scalar_bar.TitleFontSize = 8
914 scalar_bar.TitleBold = 1
915 scalar_bar.TitleItalic = 1
916 scalar_bar.TitleShadow = 1
918 scalar_bar.LabelFontFamily = 'Arial'
919 scalar_bar.LabelFontSize = 8
920 scalar_bar.LabelBold = 1
921 scalar_bar.LabelItalic = 1
922 scalar_bar.LabelShadow = 1
924 # Add the scalar bar to the view
925 pvs.GetRenderView().Representations.append(scalar_bar)
927 # Reassign the current bar
928 _current_bar = scalar_bar
934 """Get current scalar bar."""
938 def get_lookup_table(field_name, nb_components, vector_mode='Magnitude'):
939 """Get lookup table for the given field."""
940 lookup_table = pvs.GetLookupTableForArray(field_name, nb_components)
942 if vector_mode == 'Magnitude':
943 lookup_table.VectorMode = vector_mode
944 elif vector_mode == 'X':
945 lookup_table.VectorMode = 'Component'
946 lookup_table.VectorComponent = 0
947 elif vector_mode == 'Y':
948 lookup_table.VectorMode = 'Component'
949 lookup_table.VectorComponent = 1
950 elif vector_mode == 'Z':
951 lookup_table.VectorMode = 'Component'
952 lookup_table.VectorComponent = 2
954 raise ValueError("Incorrect vector mode: " + vector_mode)
956 lookup_table.Discretize = 0
957 lookup_table.ColorSpace = 'HSV'
958 if hasattr(lookup_table,"LockDataRange"):
959 lookup_table.LockDataRange = 0
960 elif hasattr(lookup_table,"LockScalarRange"):
961 lookup_table.LockScalarRange = 0
963 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
968 def get_group_mesh_name(full_group_name):
969 """Return mesh name of the group by its full name."""
970 aList = full_group_name.split('/')
972 group_name = full_group_name.split('/')[1]
976 def get_group_entity(full_group_name):
977 """Return entity type of the group by its full name."""
978 aList = full_group_name.split('/')
980 entity_name = full_group_name.split('/')[2]
981 entity = EntityType.get_type(entity_name)
985 def get_group_short_name(full_group_name):
986 """Return short name of the group by its full name."""
987 short_name = re.sub('^GRP_', '', full_group_name)
991 def get_mesh_full_names(proxy):
992 """Return all mesh names in the given proxy as a set."""
993 proxy.UpdatePipeline()
994 fields = proxy.GetProperty("FieldsTreeInfo")[::2]
995 mesh_full_names = set([item for item in fields if get_field_mesh_name(item) == get_field_short_name(item)])
996 return mesh_full_names
999 def get_group_names(extrGrps):
1000 """Return full names of all groups of the given 'ExtractGroup' filter object.
1002 group_names = filter(lambda x:x[:4]=="GRP_",list(extrGrps.GetProperty("GroupsFlagsInfo")[::2]))
1006 def get_time(proxy, timestamp_nb):
1007 """Get time value by timestamp number."""
1008 #proxy.UpdatePipeline()
1009 # Check timestamp number
1012 if (hasattr(proxy, 'TimestepValues')):
1013 timestamps = proxy.TimestepValues.GetData()
1014 elif (hasattr(proxy.Input, 'TimestepValues')):
1015 timestamps = proxy.Input.TimestepValues.GetData()
1017 length = len(timestamps)
1018 if (timestamp_nb > 0 and (timestamp_nb - 1) not in xrange(length) ) or (timestamp_nb < 0 and -timestamp_nb > length):
1019 raise ValueError("Timestamp number is out of range: " + str(timestamp_nb))
1025 if timestamp_nb > 0:
1026 return timestamps[timestamp_nb - 1]
1028 return timestamps[timestamp_nb]
1030 def create_prs(prs_type, proxy, field_entity, field_name, timestamp_nb):
1031 """Auxiliary function.
1033 Build presentation of the given type on the given field and
1035 Set the presentation properties like visu.CreatePrsForResult() do.
1038 proxy.UpdatePipeline()
1041 if prs_type == PrsTypeEnum.SCALARMAP:
1042 prs = ScalarMapOnField(proxy, field_entity, field_name, timestamp_nb)
1043 elif prs_type == PrsTypeEnum.CUTPLANES:
1044 prs = CutPlanesOnField(proxy, field_entity, field_name, timestamp_nb,
1045 orientation=Orientation.ZX)
1046 elif prs_type == PrsTypeEnum.CUTLINES:
1047 prs = CutLinesOnField(proxy, field_entity, field_name, timestamp_nb,
1048 orientation1=Orientation.XY,
1049 orientation2=Orientation.ZX)
1050 elif prs_type == PrsTypeEnum.DEFORMEDSHAPE:
1051 prs = DeformedShapeOnField(proxy, field_entity,
1052 field_name, timestamp_nb)
1053 elif prs_type == PrsTypeEnum.DEFORMEDSHAPESCALARMAP:
1054 prs = DeformedShapeAndScalarMapOnField(proxy, field_entity,
1055 field_name, timestamp_nb)
1056 elif prs_type == PrsTypeEnum.VECTORS:
1057 prs = VectorsOnField(proxy, field_entity, field_name, timestamp_nb)
1058 elif prs_type == PrsTypeEnum.PLOT3D:
1059 prs = Plot3DOnField(proxy, field_entity, field_name, timestamp_nb)
1060 elif prs_type == PrsTypeEnum.ISOSURFACES:
1061 prs = IsoSurfacesOnField(proxy, field_entity, field_name, timestamp_nb)
1062 elif prs_type == PrsTypeEnum.GAUSSPOINTS:
1063 prs = GaussPointsOnField(proxy, field_entity, field_name, timestamp_nb)
1064 elif prs_type == PrsTypeEnum.STREAMLINES:
1065 prs = StreamLinesOnField(proxy, field_entity, field_name, timestamp_nb)
1067 raise ValueError("Unexistent presentation type.")
1072 # Functions for building Post-Pro presentations
1073 def ScalarMapOnField(proxy, entity, field_name, timestamp_nb,
1074 vector_mode='Magnitude'):
1075 """Creates Scalar Map presentation on the given field.
1078 proxy: the pipeline object, containig data
1079 entity: the entity type from PrsTypeEnum
1080 field_name: the field name
1081 timestamp_nb: the number of time step (1, 2, ...)
1082 vector_mode: the mode of transformation of vector values
1083 into scalar values, applicable only if the field contains vector values.
1084 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1087 Scalar Map as representation object.
1090 proxy.UpdatePipeline()
1091 # We don't need mesh parts with no data on them
1092 if entity == EntityType.NODE:
1093 select_cells_with_data(proxy, on_points=[field_name])
1095 select_cells_with_data(proxy, on_cells=[field_name])
1098 nb_components = get_nb_components(proxy, entity, field_name)
1099 check_vector_mode(vector_mode, nb_components)
1102 time_value = get_time(proxy, timestamp_nb)
1105 pvs.GetRenderView().ViewTime = time_value
1106 pvs.UpdatePipeline(time_value, proxy)
1108 # Get Scalar Map representation object
1109 scalarmap = pvs.GetRepresentation(proxy)
1112 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
1114 # Set field range if necessary
1115 data_range = get_data_range(proxy, entity,
1116 field_name, vector_mode)
1117 if hasattr(lookup_table,"LockDataRange"):
1118 lookup_table.LockDataRange = 1
1119 elif hasattr(lookup_table,"LockScalarRange"):
1120 lookup_table.LockScalarRange = 1
1122 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1124 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1126 scalarmap.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
1127 scalarmap.LookupTable = lookup_table
1130 bar_title = field_name + ", " + str(time_value)
1131 if (nb_components > 1):
1132 bar_title += "\n" + vector_mode
1133 add_scalar_bar(field_name, nb_components, vector_mode,
1134 lookup_table, time_value)
1139 def CutPlanesOnField(proxy, entity, field_name, timestamp_nb,
1140 nb_planes=10, orientation=Orientation.YZ,
1142 displacement=0.5, vector_mode='Magnitude'):
1143 """Creates Cut Planes presentation on the given field.
1146 proxy: the pipeline object, containig data
1147 entity: the entity type from PrsTypeEnum
1148 field_name: the field name
1149 timestamp_nb: the number of time step (1, 2, ...)
1150 nb_planes: number of cutting planes
1151 orientation: cutting planes orientation in 3D space
1152 angle1: rotation of the planes in 3d space around the first axis of the
1153 selected orientation (X axis for XY, Y axis for YZ, Z axis for ZX).
1154 The angle of rotation is set in degrees. Acceptable range: [-45, 45].
1155 angle2: rotation of the planes in 3d space around the second axis of the
1156 selected orientation. Acceptable range: [-45, 45].
1157 displacement: the displacement of the planes into one or another side
1158 vector_mode: the mode of transformation of vector values
1159 into scalar values, applicable only if the field contains vector values.
1160 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1163 Cut Planes as representation object.
1166 proxy.UpdatePipeline()
1167 if entity == EntityType.NODE:
1168 select_cells_with_data(proxy, on_points=[field_name])
1170 select_cells_with_data(proxy, on_cells=[field_name])
1173 nb_components = get_nb_components(proxy, entity, field_name)
1174 check_vector_mode(vector_mode, nb_components)
1177 time_value = get_time(proxy, timestamp_nb)
1180 pvs.GetRenderView().ViewTime = time_value
1181 pvs.UpdatePipeline(time_value, proxy)
1183 # Create slice filter
1184 slice_filter = pvs.Slice(proxy)
1185 slice_filter.SliceType = "Plane"
1187 # Set cut planes normal
1188 normal = get_normal_by_orientation(orientation,
1189 radians(angle1), radians(angle2))
1190 slice_filter.SliceType.Normal = normal
1192 # Set cut planes positions
1193 positions = get_positions(nb_planes, normal,
1194 get_bounds(proxy), displacement)
1195 slice_filter.SliceOffsetValues = positions
1197 # Get Cut Planes representation object
1198 cut_planes = pvs.GetRepresentation(slice_filter)
1201 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
1203 # Set field range if necessary
1204 data_range = get_data_range(proxy, entity,
1205 field_name, vector_mode)
1207 if hasattr(lookup_table,"LockDataRange"):
1208 lookup_table.LockDataRange = 1
1209 elif hasattr(lookup_table,"LockScalarRange"):
1210 lookup_table.LockScalarRange = 1
1212 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1214 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1217 cut_planes.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
1218 cut_planes.LookupTable = lookup_table
1221 add_scalar_bar(field_name, nb_components,
1222 vector_mode, lookup_table, time_value)
1227 def CutLinesOnField(proxy, entity, field_name, timestamp_nb,
1229 orientation1=Orientation.XY,
1230 base_angle1=0, base_angle2=0,
1231 orientation2=Orientation.YZ,
1232 cut_angle1=0, cut_angle2=0,
1233 displacement1=0.5, displacement2=0.5,
1234 generate_curves=False,
1235 vector_mode='Magnitude'):
1236 """Creates Cut Lines presentation on the given field.
1239 proxy: the pipeline object, containig data
1240 entity: the entity type from PrsTypeEnum
1241 field_name: the field name
1242 timestamp_nb: the number of time step (1, 2, ...)
1243 nb_lines: number of lines
1244 orientation1: base plane orientation in 3D space
1245 base_angle1: rotation of the base plane in 3d space around the first
1246 axis of the orientation1 (X axis for XY, Y axis for YZ, Z axis for ZX).
1247 The angle of rotation is set in degrees. Acceptable range: [-45, 45].
1248 base_angle2: rotation of the base plane in 3d space around the second
1249 axis of the orientation1. Acceptable range: [-45, 45].
1250 orientation2: cutting planes orientation in 3D space
1251 cut_angle1: rotation of the cut planes in 3d space around the first
1252 axis of the orientation2. Acceptable range: [-45, 45].
1253 cut_angle2: rotation of the cuting planes in 3d space around the second
1254 axis of the orientation2. Acceptable range: [-45, 45].
1255 displacement1: base plane displacement
1256 displacement2: cutting planes displacement
1257 generate_curves: if true, 'PlotOverLine' filter will be created
1259 vector_mode: the mode of transformation of vector values
1260 into scalar values, applicable only if the field contains vector values.
1261 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1264 Cut Lines as representation object if generate_curves == False,
1265 (Cut Lines as representation object, list of 'PlotOverLine') otherwise
1268 proxy.UpdatePipeline()
1269 if entity == EntityType.NODE:
1270 select_cells_with_data(proxy, on_points=[field_name])
1272 select_cells_with_data(proxy, on_cells=[field_name])
1275 nb_components = get_nb_components(proxy, entity, field_name)
1276 check_vector_mode(vector_mode, nb_components)
1279 time_value = get_time(proxy, timestamp_nb)
1282 pvs.GetRenderView().ViewTime = time_value
1283 pvs.UpdatePipeline(time_value, proxy)
1286 base_plane = pvs.Slice(proxy)
1287 base_plane.SliceType = "Plane"
1289 # Set base plane normal
1290 base_normal = get_normal_by_orientation(orientation1,
1291 radians(base_angle1),
1292 radians(base_angle2))
1293 base_plane.SliceType.Normal = base_normal
1295 # Set base plane position
1296 base_position = get_positions(1, base_normal,
1297 get_bounds(proxy), displacement1)
1298 base_plane.SliceOffsetValues = base_position
1301 base_plane.UpdatePipeline()
1302 if (base_plane.GetDataInformation().GetNumberOfCells() == 0):
1305 # Create cutting planes
1306 cut_planes = pvs.Slice(base_plane)
1307 cut_planes.SliceType = "Plane"
1309 # Set cutting planes normal and get positions
1310 cut_normal = get_normal_by_orientation(orientation2,
1311 radians(cut_angle1),
1312 radians(cut_angle2))
1313 cut_planes.SliceType.Normal = cut_normal
1315 # Set cutting planes position
1316 cut_positions = get_positions(nb_lines, cut_normal,
1317 get_bounds(base_plane), displacement2)
1323 for pos in cut_positions:
1324 # Get points for plot over line objects
1325 cut_planes.SliceOffsetValues = pos
1326 cut_planes.UpdatePipeline()
1327 bounds = get_bounds(cut_planes)
1328 point1 = [bounds[0], bounds[2], bounds[4]]
1329 point2 = [bounds[1], bounds[3], bounds[5]]
1331 # Create plot over line filter
1332 pol = pvs.PlotOverLine(cut_planes,
1333 Source="High Resolution Line Source")
1334 pvs.RenameSource('Y' + str(index), pol)
1335 pol.Source.Point1 = point1
1336 pol.Source.Point2 = point2
1337 pol.UpdatePipeline()
1342 cut_planes.SliceOffsetValues = cut_positions
1343 cut_planes.UpdatePipeline()
1345 # Get Cut Lines representation object
1346 cut_lines = pvs.GetRepresentation(cut_planes)
1349 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
1351 # Set field range if necessary
1352 data_range = get_data_range(proxy, entity,
1353 field_name, vector_mode)
1354 if hasattr(lookup_table,"LockDataRange"):
1355 lookup_table.LockDataRange = 1
1356 elif hasattr(lookup_table,"LockScalarRange"):
1357 lookup_table.LockScalarRange = 1
1359 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1361 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1364 cut_lines.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
1365 cut_lines.LookupTable = lookup_table
1367 # Set wireframe represenatation mode
1368 cut_lines.Representation = 'Wireframe'
1371 add_scalar_bar(field_name, nb_components,
1372 vector_mode, lookup_table, time_value)
1375 # If curves were generated return tuple (cut lines, list of curves)
1377 result = cut_lines, curves
1382 def CutSegmentOnField(proxy, entity, field_name, timestamp_nb,
1383 point1, point2, vector_mode='Magnitude'):
1384 """Creates Cut Segment presentation on the given field.
1387 proxy: the pipeline object, containig data
1388 entity: the entity type from PrsTypeEnum
1389 field_name: the field name
1390 timestamp_nb: the number of time step (1, 2, ...)
1391 point1: set the first point of the segment (as [x, y, z])
1392 point1: set the second point of the segment (as [x, y, z])
1393 vector_mode: the mode of transformation of vector values
1394 into scalar values, applicable only if the field contains vector values.
1395 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1398 Cut Segment as 3D representation object.
1401 proxy.UpdatePipeline()
1402 if entity == EntityType.NODE:
1403 select_cells_with_data(proxy, on_points=[field_name])
1405 select_cells_with_data(proxy, on_cells=[field_name])
1408 nb_components = get_nb_components(proxy, entity, field_name)
1409 check_vector_mode(vector_mode, nb_components)
1412 time_value = get_time(proxy, timestamp_nb)
1415 pvs.GetRenderView().ViewTime = time_value
1416 pvs.UpdatePipeline(time_value, proxy)
1418 # Create plot over line filter
1419 pol = pvs.PlotOverLine(proxy, Source="High Resolution Line Source")
1420 pol.Source.Point1 = point1
1421 pol.Source.Point2 = point2
1422 pol.UpdatePipeline()
1424 # Get Cut Segment representation object
1425 cut_segment = pvs.GetRepresentation(pol)
1428 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
1430 # Set field range if necessary
1431 data_range = get_data_range(proxy, entity,
1432 field_name, vector_mode)
1433 if hasattr(lookup_table,"LockDataRange"):
1434 lookup_table.LockDataRange = 1
1435 elif hasattr(lookup_table,"LockScalarRange"):
1436 lookup_table.LockScalarRange = 1
1438 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1440 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1443 cut_segment.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
1444 cut_segment.LookupTable = lookup_table
1446 # Set wireframe represenatation mode
1447 cut_segment.Representation = 'Wireframe'
1450 add_scalar_bar(field_name, nb_components,
1451 vector_mode, lookup_table, time_value)
1456 def VectorsOnField(proxy, entity, field_name, timestamp_nb,
1458 glyph_pos=GlyphPos.TAIL, glyph_type='2D Glyph',
1459 is_colored=False, vector_mode='Magnitude'):
1460 """Creates Vectors presentation on the given field.
1463 proxy: the pipeline object, containig data
1464 entity: the entity type from PrsTypeEnum
1465 field_name: the field name
1466 timestamp_nb: the number of time step (1, 2, ...)
1467 scale_factor: scale factor
1468 glyph_pos: the position of glyphs
1469 glyph_type: the type of glyphs
1470 is_colored: this option allows to color the presentation according to
1471 the corresponding data array values
1472 vector_mode: the mode of transformation of vector values
1473 into scalar values, applicable only if the field contains vector values.
1474 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1477 Vectors as representation object.
1480 proxy.UpdatePipeline()
1481 if entity == EntityType.NODE:
1482 select_cells_with_data(proxy, on_points=[field_name])
1484 select_cells_with_data(proxy, on_cells=[field_name])
1487 nb_components = get_nb_components(proxy, entity, field_name)
1488 check_vector_mode(vector_mode, nb_components)
1491 time_value = get_time(proxy, timestamp_nb)
1494 pvs.GetRenderView().ViewTime = time_value
1495 pvs.UpdatePipeline(time_value, proxy)
1497 # Extract only groups with data for the field
1501 if is_data_on_cells(proxy, field_name):
1502 cell_centers = pvs.CellCenters(source)
1503 cell_centers.VertexCells = 1
1504 source = cell_centers
1506 vector_array = field_name
1507 # If the given vector array has only 2 components, add the third one
1508 if nb_components == 2:
1509 calc = get_add_component_calc(source, EntityType.NODE, field_name)
1510 vector_array = calc.ResultArrayName
1514 glyph = pvs.Glyph(source)
1515 glyph.Vectors = vector_array
1516 glyph.ScaleMode = 'vector'
1517 #glyph.MaskPoints = 0
1520 glyph.GlyphType = glyph_type
1521 if glyph_type == '2D Glyph':
1522 glyph.GlyphType.GlyphType = 'Arrow'
1523 elif glyph_type == 'Cone':
1524 glyph.GlyphType.Resolution = 7
1525 glyph.GlyphType.Height = 2
1526 glyph.GlyphType.Radius = 0.2
1528 # Set glyph position if possible
1529 if glyph.GlyphType.GetProperty("Center"):
1530 if (glyph_pos == GlyphPos.TAIL):
1531 glyph.GlyphType.Center = [0.5, 0.0, 0.0]
1532 elif (glyph_pos == GlyphPos.HEAD):
1533 glyph.GlyphType.Center = [-0.5, 0.0, 0.0]
1534 elif (glyph_pos == GlyphPos.CENTER):
1535 glyph.GlyphType.Center = [0.0, 0.0, 0.0]
1537 if scale_factor is not None:
1538 glyph.ScaleFactor = scale_factor
1540 def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE,
1541 proxy, entity, field_name)
1542 glyph.ScaleFactor = def_scale
1544 glyph.UpdatePipeline()
1546 # Get Vectors representation object
1547 vectors = pvs.GetRepresentation(glyph)
1550 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
1552 # Set field range if necessary
1553 data_range = get_data_range(proxy, entity,
1554 field_name, vector_mode)
1555 if hasattr(lookup_table,"LockDataRange"):
1556 lookup_table.LockDataRange = 1
1557 elif hasattr(lookup_table,"LockScalarRange"):
1558 lookup_table.LockScalarRange = 1
1560 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1562 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1566 vectors.ColorArrayName = (EntityType.get_pvtype(entity), 'GlyphVector')
1568 vectors.ColorArrayName = (None, '')
1569 vectors.LookupTable = lookup_table
1571 vectors.LineWidth = 1.0
1573 # Set wireframe represenatation mode
1574 vectors.Representation = 'Wireframe'
1577 add_scalar_bar(field_name, nb_components,
1578 vector_mode, lookup_table, time_value)
1583 def DeformedShapeOnField(proxy, entity, field_name,
1585 scale_factor=None, is_colored=False,
1586 vector_mode='Magnitude'):
1587 """Creates Defromed Shape presentation on the given field.
1590 proxy: the pipeline object, containig data
1591 entity: the entity type from PrsTypeEnum
1592 field_name: the field name
1593 timestamp_nb: the number of time step (1, 2, ...)
1594 scale_factor: scale factor of the deformation
1595 is_colored: this option allows to color the presentation according to
1596 the corresponding data array values
1597 vector_mode: the mode of transformation of vector values
1598 into scalar values, applicable only if the field contains vector values.
1599 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1602 Defromed Shape as representation object.
1605 proxy.UpdatePipeline()
1606 # We don't need mesh parts with no data on them
1607 if entity == EntityType.NODE:
1608 select_cells_with_data(proxy, on_points=[field_name])
1610 select_cells_with_data(proxy, on_cells=[field_name])
1613 nb_components = get_nb_components(proxy, entity, field_name)
1614 check_vector_mode(vector_mode, nb_components)
1617 time_value = get_time(proxy, timestamp_nb)
1620 pvs.GetRenderView().ViewTime = time_value
1621 pvs.UpdatePipeline(time_value, proxy)
1624 source = pvs.MergeBlocks(proxy)
1625 pvs.UpdatePipeline()
1627 # Cell data to point data
1628 if is_data_on_cells(proxy, field_name):
1629 cell_to_point = pvs.CellDatatoPointData()
1630 cell_to_point.PassCellData = 1
1631 source = cell_to_point
1633 vector_array = field_name
1634 # If the given vector array has only 2 components, add the third one
1635 if nb_components == 2:
1636 calc = get_add_component_calc(source, EntityType.NODE, field_name)
1637 vector_array = calc.ResultArrayName
1641 warp_vector = pvs.WarpByVector(source)
1642 warp_vector.Vectors = [vector_array]
1643 if scale_factor is not None:
1644 warp_vector.ScaleFactor = scale_factor
1646 def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE,
1647 proxy, entity, field_name)
1648 warp_vector.ScaleFactor = def_scale
1650 # Get Deformed Shape representation object
1651 defshape = pvs.GetRepresentation(warp_vector)
1654 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
1656 # Set field range if necessary
1657 data_range = get_data_range(proxy, entity,
1658 field_name, vector_mode)
1659 if hasattr(lookup_table,"LockDataRange"):
1660 lookup_table.LockDataRange = 1
1661 elif hasattr(lookup_table,"LockScalarRange"):
1662 lookup_table.LockScalarRange = 1
1664 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1666 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1670 defshape.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
1672 defshape.ColorArrayName = (None, '')
1673 defshape.LookupTable = lookup_table
1675 # Set wireframe represenatation mode
1676 defshape.Representation = 'Wireframe'
1679 add_scalar_bar(field_name, nb_components,
1680 vector_mode, lookup_table, time_value)
1685 def DeformedShapeAndScalarMapOnField(proxy, entity, field_name,
1689 scalar_field_name=None,
1690 vector_mode='Magnitude'):
1691 """Creates Defromed Shape And Scalar Map presentation on the given field.
1694 proxy: the pipeline object, containig data
1695 entity: the entity type from PrsTypeEnum
1696 field_name: the field name
1697 timestamp_nb: the number of time step (1, 2, ...)
1698 scale_factor: scale factor of the deformation
1699 scalar_entity: scalar field entity
1700 scalar_field_name: scalar field, i.e. the field for coloring
1701 vector_mode: the mode of transformation of vector values
1702 into scalar values, applicable only if the field contains vector values.
1703 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1706 Defromed Shape And Scalar Map as representation object.
1709 proxy.UpdatePipeline()
1710 # We don't need mesh parts with no data on them
1714 if entity == EntityType.NODE:
1715 on_points.append(field_name)
1717 on_cells.append(field_name)
1719 if scalar_entity and scalar_field_name:
1720 if scalar_entity == EntityType.NODE:
1721 on_points.append(scalar_field_name)
1723 on_cells.append(scalar_field_name)
1725 nb_components = get_nb_components(proxy, entity, field_name)
1728 select_cells_with_data(proxy, on_points, on_cells)
1731 check_vector_mode(vector_mode, nb_components)
1734 time_value = get_time(proxy, timestamp_nb)
1737 pvs.GetRenderView().ViewTime = time_value
1738 pvs.UpdatePipeline(time_value, proxy)
1740 # Set scalar field by default
1741 scalar_field_entity = scalar_entity
1742 scalar_field = scalar_field_name
1743 if (scalar_field_entity is None) or (scalar_field is None):
1744 scalar_field_entity = entity
1745 scalar_field = field_name
1748 source = pvs.MergeBlocks(proxy)
1749 pvs.UpdatePipeline()
1751 # Cell data to point data
1752 if is_data_on_cells(proxy, field_name):
1753 cell_to_point = pvs.CellDatatoPointData(source)
1754 cell_to_point.PassCellData = 1
1755 source = cell_to_point
1757 vector_array = field_name
1758 # If the given vector array has only 2 components, add the third one
1759 if nb_components == 2:
1760 calc = get_add_component_calc(source, EntityType.NODE, field_name)
1761 vector_array = calc.ResultArrayName
1765 warp_vector = pvs.WarpByVector(source)
1766 warp_vector.Vectors = [vector_array]
1767 if scale_factor is not None:
1768 warp_vector.ScaleFactor = scale_factor
1770 def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE,
1771 proxy, entity, field_name)
1772 warp_vector.ScaleFactor = def_scale
1774 # Get Defromed Shape And Scalar Map representation object
1775 defshapemap = pvs.GetRepresentation(warp_vector)
1778 lookup_table = get_lookup_table(scalar_field, nb_components, vector_mode)
1780 # Set field range if necessary
1781 data_range = get_data_range(proxy, scalar_field_entity,
1782 scalar_field, vector_mode)
1783 if hasattr(lookup_table,"LockDataRange"):
1784 lookup_table.LockDataRange = 1
1785 elif hasattr(lookup_table,"LockScalarRange"):
1786 lookup_table.LockScalarRange = 1
1788 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1790 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1793 defshapemap.ColorArrayName = (EntityType.get_pvtype(scalar_field_entity), scalar_field)
1794 defshapemap.LookupTable = lookup_table
1797 add_scalar_bar(field_name, nb_components,
1798 vector_mode, lookup_table, time_value)
1803 def Plot3DOnField(proxy, entity, field_name, timestamp_nb,
1804 orientation=Orientation.AUTO,
1806 position=0.5, is_relative=True,
1808 is_contour=False, nb_contours=32,
1809 vector_mode='Magnitude'):
1810 """Creates Plot 3D presentation on the given field.
1813 proxy: the pipeline object, containig data
1814 entity: the entity type from PrsTypeEnum
1815 field_name: the field name
1816 timestamp_nb: the number of time step (1, 2, ...)
1817 orientation: the cut plane plane orientation in 3D space, if
1818 the input is planar - will not be taken into account
1819 angle1: rotation of the cut plane in 3d space around the first axis
1820 of the selected orientation (X axis for XY, Y axis for YZ,
1822 The angle of rotation is set in degrees. Acceptable range: [-45, 45].
1823 angle2: rotation of the cut plane in 3d space around the second axis
1824 of the selected orientation. Acceptable range: [-45, 45].
1825 position: position of the cut plane in the object (ranging from 0 to 1).
1826 The value 0.5 corresponds to cutting by halves.
1827 is_relative: defines if the cut plane position is relative or absolute
1828 scale_factor: deformation scale factor
1829 is_contour: if True - Plot 3D will be represented with a set of contours,
1830 otherwise - Plot 3D will be represented with a smooth surface
1831 nb_contours: number of contours, applied if is_contour is True
1832 vector_mode: the mode of transformation of vector values
1833 into scalar values, applicable only if the field contains vector values.
1834 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
1837 Plot 3D as representation object.
1840 proxy.UpdatePipeline()
1841 # We don't need mesh parts with no data on them
1842 if entity == EntityType.NODE:
1843 select_cells_with_data(proxy, on_points=[field_name])
1845 select_cells_with_data(proxy, on_cells=[field_name])
1848 nb_components = get_nb_components(proxy, entity, field_name)
1849 check_vector_mode(vector_mode, nb_components)
1852 time_value = get_time(proxy, timestamp_nb)
1855 pvs.GetRenderView().ViewTime = time_value
1856 pvs.UpdatePipeline(time_value, proxy)
1859 merge_blocks = pvs.MergeBlocks(proxy)
1860 merge_blocks.UpdatePipeline()
1866 # Define orientation if necessary (auto mode)
1867 plane_orientation = orientation
1868 if (orientation == Orientation.AUTO):
1869 plane_orientation = get_orientation(proxy)
1871 # Get cutting plane normal
1874 if (not is_planar_input(proxy)):
1875 normal = get_normal_by_orientation(plane_orientation,
1876 radians(angle1), radians(angle2))
1878 # Create slice filter
1879 slice_filter = pvs.Slice(merge_blocks)
1880 slice_filter.SliceType = "Plane"
1882 # Set cutting plane normal
1883 slice_filter.SliceType.Normal = normal
1885 # Set cutting plane position
1887 base_position = get_positions(1, normal,
1888 get_bounds(proxy), position)
1889 slice_filter.SliceOffsetValues = base_position
1891 slice_filter.SliceOffsetValues = position
1893 slice_filter.UpdatePipeline()
1894 poly_data = slice_filter
1896 normal = get_normal_by_orientation(plane_orientation, 0, 0)
1900 if not poly_data or poly_data.GetDataInformation().GetNumberOfCells() == 0:
1901 geometry_filter = pvs.GeometryFilter(merge_blocks)
1902 poly_data = geometry_filter
1903 use_normal = 1 # TODO(MZN): workaround
1909 if is_data_on_cells(poly_data, field_name):
1910 # Cell data to point data
1911 cell_to_point = pvs.CellDatatoPointData(poly_data)
1912 cell_to_point.PassCellData = 1
1913 source = cell_to_point
1915 scalars = ['POINTS', field_name]
1917 # Transform vector array to scalar array if necessary
1918 if (nb_components > 1):
1919 calc = get_calc_magnitude(source, EntityType.NODE, field_name)
1920 scalars = ['POINTS', calc.ResultArrayName]
1924 warp_scalar = pvs.WarpByScalar(source)
1925 warp_scalar.Scalars = scalars
1926 warp_scalar.Normal = normal
1927 warp_scalar.UseNormal = use_normal
1928 if scale_factor is not None:
1929 warp_scalar.ScaleFactor = scale_factor
1931 def_scale = get_default_scale(PrsTypeEnum.PLOT3D,
1932 proxy, entity, field_name)
1933 warp_scalar.ScaleFactor = def_scale
1935 warp_scalar.UpdatePipeline()
1936 source = warp_scalar
1940 contour = pvs.Contour(warp_scalar)
1941 contour.PointMergeMethod = "Uniform Binning"
1942 contour.ContourBy = ['POINTS', field_name]
1943 scalar_range = get_data_range(proxy, entity,
1944 field_name, vector_mode)
1945 contour.Isosurfaces = get_contours(scalar_range, nb_contours)
1946 contour.UpdatePipeline()
1949 # Get Plot 3D representation object
1950 plot3d = pvs.GetRepresentation(source)
1953 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
1955 # Set field range if necessary
1956 data_range = get_data_range(proxy, entity,
1957 field_name, vector_mode)
1958 if hasattr(lookup_table,"LockDataRange"):
1959 lookup_table.LockDataRange = 1
1960 elif hasattr(lookup_table,"LockScalarRange"):
1961 lookup_table.LockScalarRange = 1
1963 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
1965 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
1968 plot3d.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
1969 plot3d.LookupTable = lookup_table
1972 add_scalar_bar(field_name, nb_components,
1973 vector_mode, lookup_table, time_value)
1978 def IsoSurfacesOnField(proxy, entity, field_name, timestamp_nb,
1979 custom_range=None, nb_surfaces=10,
1980 is_colored=True, color=None, vector_mode='Magnitude'):
1981 """Creates Iso Surfaces presentation on the given field.
1984 proxy: the pipeline object, containig data
1985 entity: the entity type from PrsTypeEnum
1986 field_name: the field name
1987 timestamp_nb: the number of time step (1, 2, ...)
1988 custom_range: scalar range, if undefined the source range will be applied
1989 nb_surfaces: number of surfaces, which will be generated
1990 is_colored: this option allows to color the presentation according to
1991 the corresponding data array values. If False - the presentation will
1993 color: defines the presentation color as [R, G, B] triple. Taken into
1994 account only if is_colored is False.
1995 vector_mode: the mode of transformation of vector values
1996 into scalar values, applicable only if the field contains vector values.
1997 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
2000 Iso Surfaces as representation object.
2003 proxy.UpdatePipeline()
2004 # We don't need mesh parts with no data on them
2005 if entity == EntityType.NODE:
2006 select_cells_with_data(proxy, on_points=[field_name])
2008 select_cells_with_data(proxy, on_cells=[field_name])
2011 nb_components = get_nb_components(proxy, entity, field_name)
2012 check_vector_mode(vector_mode, nb_components)
2015 time_value = get_time(proxy, timestamp_nb)
2018 pvs.GetRenderView().ViewTime = time_value
2019 pvs.UpdatePipeline(time_value, proxy)
2022 source = pvs.MergeBlocks(proxy)
2023 pvs.UpdatePipeline()
2025 # Transform cell data into point data if necessary
2026 if is_data_on_cells(proxy, field_name):
2027 cell_to_point = pvs.CellDatatoPointData(source)
2028 cell_to_point.PassCellData = 1
2029 source = cell_to_point
2031 contour_by = ['POINTS', field_name]
2033 # Transform vector array to scalar array if necessary
2034 if (nb_components > 1):
2035 calc = get_calc_magnitude(source, EntityType.NODE, field_name)
2036 contour_by = ['POINTS', calc.ResultArrayName]
2039 # Contour filter settings
2040 contour = pvs.Contour(source)
2041 contour.ComputeScalars = 1
2042 contour.ContourBy = contour_by
2045 scalar_range = custom_range
2046 if (scalar_range is None):
2047 scalar_range = get_data_range(proxy, entity,
2048 field_name, cut_off=True)
2050 # Get contour values for the range
2051 surfaces = get_contours(scalar_range, nb_surfaces)
2053 # Set contour values
2054 contour.Isosurfaces = surfaces
2056 # Get Iso Surfaces representation object
2057 isosurfaces = pvs.GetRepresentation(contour)
2060 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
2062 # Set field range if necessary
2063 data_range = get_data_range(proxy, entity,
2064 field_name, vector_mode)
2065 if hasattr(lookup_table,"LockDataRange"):
2066 lookup_table.LockDataRange = 1
2067 elif hasattr(lookup_table,"LockScalarRange"):
2068 lookup_table.LockScalarRange = 1
2070 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
2072 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
2074 # Set display properties
2076 isosurfaces.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
2078 isosurfaces.ColorArrayName = (None, '')
2080 isosurfaces.DiffuseColor = color
2081 isosurfaces.LookupTable = lookup_table
2084 add_scalar_bar(field_name, nb_components,
2085 vector_mode, lookup_table, time_value)
2090 def GaussPointsOnField(proxy, entity, field_name,
2092 is_deformed=True, scale_factor=None,
2093 is_colored=True, color=None,
2094 primitive=GaussType.SPRITE,
2095 is_proportional=True,
2097 multiplier=None, vector_mode='Magnitude'):
2098 """Creates Gauss Points on the given field.
2102 proxy: the pipeline object, containig data
2103 entity: the field entity type from PrsTypeEnum
2104 field_name: the field name
2105 timestamp_nb: the number of time step (1, 2, ...)
2106 is_deformed: defines whether the Gauss Points will be deformed or not
2107 scale_factor -- the scale factor for deformation. Will be taken into
2108 account only if is_deformed is True.
2109 If not passed by user, default scale will be computed.
2110 is_colored -- defines whether the Gauss Points will be multicolored,
2111 using the corresponding data values
2112 color: defines the presentation color as [R, G, B] triple. Taken into
2113 account only if is_colored is False.
2114 primitive: primitive type from GaussType
2115 is_proportional: if True, the size of primitives will depends on
2116 the gauss point value
2117 max_pixel_size: the maximum sizr of the Gauss Points primitive in pixels
2118 multiplier: coefficient between data values and the size of primitives
2119 If not passed by user, default scale will be computed.
2120 vector_mode: the mode of transformation of vector values into
2121 scalar values, applicable only if the field contains vector values.
2122 Possible modes: 'Magnitude' - vector module;
2123 'X', 'Y', 'Z' - vector components.
2126 Gauss Points as representation object.
2129 proxy.UpdatePipeline()
2130 # We don't need mesh parts with no data on them
2131 on_gauss = select_cells_with_data(proxy, on_gauss=[field_name])
2133 if entity == EntityType.NODE:
2134 select_cells_with_data(proxy, on_points=[field_name])
2136 select_cells_with_data(proxy, on_cells=[field_name])
2139 nb_components = get_nb_components(proxy, entity, field_name)
2140 check_vector_mode(vector_mode, nb_components)
2143 time_value = get_time(proxy, timestamp_nb)
2146 pvs.GetRenderView().ViewTime = time_value
2147 pvs.UpdatePipeline(time_value, proxy)
2151 # If no quadrature point array is passed, use cell centers
2153 generate_qp = pvs.GenerateQuadraturePoints(source)
2154 generate_qp.QuadratureSchemeDef = ['CELLS', 'ELGA@0']
2155 source = generate_qp
2158 cell_centers = pvs.CellCenters(source)
2159 cell_centers.VertexCells = 1
2160 source = cell_centers
2162 source.UpdatePipeline()
2164 # Check if deformation enabled
2165 if is_deformed and nb_components > 1:
2166 vector_array = field_name
2167 # If the given vector array has only 2 components, add the third one
2168 if nb_components == 2:
2169 calc = get_add_component_calc(source, EntityType.NODE, field_name)
2170 vector_array = calc.ResultArrayName
2174 warp_vector = pvs.WarpByVector(source)
2175 warp_vector.Vectors = [vector_array]
2176 if scale_factor is not None:
2177 warp_vector.ScaleFactor = scale_factor
2179 def_scale = get_default_scale(PrsTypeEnum.DEFORMEDSHAPE, proxy,
2181 warp_vector.ScaleFactor = def_scale
2182 warp_vector.UpdatePipeline()
2183 source = warp_vector
2185 # Get Gauss Points representation object
2186 gausspnt = pvs.GetRepresentation(source)
2189 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
2191 # Set field range if necessary
2192 data_range = get_data_range(proxy, entity,
2193 field_name, vector_mode)
2194 if hasattr(lookup_table,"LockDataRange"):
2195 lookup_table.LockDataRange = 1
2196 elif hasattr(lookup_table,"LockScalarRange"):
2197 lookup_table.LockScalarRange = 1
2199 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
2201 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
2203 # Set display properties
2205 gausspnt.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
2207 gausspnt.ColorArrayName = (None, '')
2209 gausspnt.DiffuseColor = color
2211 gausspnt.LookupTable = lookup_table
2214 add_scalar_bar(field_name, nb_components,
2215 vector_mode, lookup_table, time_value)
2217 # Set point sprite representation
2218 gausspnt.Representation = 'Point Sprite'
2220 # Point sprite settings
2221 gausspnt.InterpolateScalarsBeforeMapping = 0
2222 gausspnt.MaxPixelSize = max_pixel_size
2225 gausspnt.RenderMode = GaussType.get_mode(primitive)
2227 #if primitive == GaussType.SPRITE:
2229 # TODO(MZN): replace with pvsimple high-level interface
2230 # texture = sm.CreateProxy("textures", "SpriteTexture")
2231 # alphamprop = texture.GetProperty("AlphaMethod")
2232 # alphamprop.SetElement(0, 2) # Clamp
2233 # alphatprop = texture.GetProperty("AlphaThreshold")
2234 # alphatprop.SetElement(0, 63)
2235 # maxprop = texture.GetProperty("Maximum")
2236 # maxprop.SetElement(0, 255)
2237 # texture.UpdateVTKObjects()
2239 # gausspnt.Texture = texture
2240 #gausspnt.Texture.AlphaMethod = 'Clamp'
2241 #gausspnt.Texture.AlphaThreshold = 63
2242 #gausspnt.Texture.Maximum= 255
2244 # Proportional radius
2245 gausspnt.RadiusUseScalarRange = 0
2246 gausspnt.RadiusIsProportional = 0
2250 if mult is None and data_range[1] != 0:
2251 mult = abs(0.1 / data_range[1])
2253 gausspnt.RadiusScalarRange = data_range
2254 gausspnt.RadiusTransferFunctionEnabled = 1
2255 gausspnt.RadiusMode = 'Scalar'
2256 gausspnt.RadiusArray = ['POINTS', field_name]
2257 if nb_components > 1:
2258 v_comp = get_vector_component(vector_mode)
2259 gausspnt.RadiusVectorComponent = v_comp
2260 gausspnt.RadiusTransferFunctionMode = 'Table'
2261 gausspnt.RadiusScalarRange = data_range
2262 gausspnt.RadiusUseScalarRange = 1
2263 if mult is not None:
2264 gausspnt.RadiusIsProportional = 1
2265 gausspnt.RadiusProportionalFactor = mult
2267 gausspnt.RadiusTransferFunctionEnabled = 0
2268 gausspnt.RadiusMode = 'Constant'
2269 gausspnt.RadiusArray = ['POINTS', 'Constant Radius']
2273 def GaussPointsOnField1(proxy, entity, field_name,
2275 is_colored=True, color=None,
2276 primitive=GaussType.SPHERE,
2277 is_proportional=True,
2280 vector_mode='Magnitude'):
2281 """Creates Gauss Points on the given field. Use GaussPoints() Paraview interface.
2284 proxy: the pipeline object, containig data
2285 entity: the field entity type from PrsTypeEnum
2286 field_name: the field name
2287 timestamp_nb: the number of time step (1, 2, ...)
2288 is_colored -- defines whether the Gauss Points will be multicolored,
2289 using the corresponding data values
2290 color: defines the presentation color as [R, G, B] triple. Taken into
2291 account only if is_colored is False.
2292 primitive: primitive type from GaussType
2293 is_proportional: if True, the size of primitives will depends on
2294 the gauss point value
2295 max_pixel_size: the maximum sizr of the Gauss Points primitive in pixels
2296 multiplier: coefficient between data values and the size of primitives
2297 If not passed by user, default scale will be computed.
2298 vector_mode: the mode of transformation of vector values into
2299 scalar values, applicable only if the field contains vector values.
2300 Possible modes: 'Magnitude' - vector module;
2301 'X', 'Y', 'Z' - vector components.
2304 Gauss Points as representation object.
2307 proxy.UpdatePipeline()
2308 select_cells_with_data(proxy, on_gauss=[field_name])
2310 nb_components = get_nb_components(proxy, entity, field_name)
2313 time_value = get_time(proxy, timestamp_nb)
2316 pvs.GetRenderView().ViewTime = time_value
2317 proxy.UpdatePipeline(time=time_value)
2319 # Create Gauss Points object
2320 source = pvs.GaussPoints(proxy)
2321 source.UpdatePipeline()
2323 # Get Gauss Points representation object
2324 gausspnt = pvs.GetRepresentation(source)
2327 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
2329 # Set field range if necessary
2330 data_range = get_data_range(proxy, entity,
2331 field_name, vector_mode)
2332 if hasattr(lookup_table,"LockDataRange"):
2333 lookup_table.LockDataRange = 1
2334 elif hasattr(lookup_table,"LockScalarRange"):
2335 lookup_table.LockScalarRange = 1
2337 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
2339 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
2341 # Set display properties
2343 gausspnt.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
2345 gausspnt.ColorArrayName = (None, '')
2347 gausspnt.DiffuseColor = color
2349 gausspnt.LookupTable = lookup_table
2352 add_scalar_bar(field_name, nb_components,
2353 vector_mode, lookup_table, time_value)
2355 # Set point sprite representation
2356 gausspnt.Representation = 'Point Sprite'
2358 # Point sprite settings
2359 gausspnt.InterpolateScalarsBeforeMapping = 0
2360 gausspnt.MaxPixelSize = max_pixel_size
2363 gausspnt.RenderMode = GaussType.get_mode(primitive)
2365 #if primitive == GaussType.SPRITE:
2367 # TODO(MZN): replace with pvsimple high-level interface
2368 # texture = sm.CreateProxy("textures", "SpriteTexture")
2369 # alphamprop = texture.GetProperty("AlphaMethod")
2370 # alphamprop.SetElement(0, 2) # Clamp
2371 # alphatprop = texture.GetProperty("AlphaThreshold")
2372 # alphatprop.SetElement(0, 63)
2373 # maxprop = texture.GetProperty("Maximum")
2374 # maxprop.SetElement(0, 255)
2375 # texture.UpdateVTKObjects()
2377 # gausspnt.Texture = texture
2378 #gausspnt.Texture.AlphaMethod = 'Clamp'
2379 #gausspnt.Texture.AlphaThreshold = 63
2380 #gausspnt.Texture.Maximum= 255
2382 # Proportional radius
2383 gausspnt.RadiusUseScalarRange = 0
2384 gausspnt.RadiusIsProportional = 0
2388 if mult is None and data_range[1] != 0:
2389 mult = abs(0.1 / data_range[1])
2391 gausspnt.RadiusScalarRange = data_range
2392 gausspnt.RadiusTransferFunctionEnabled = 1
2393 gausspnt.RadiusMode = 'Scalar'
2394 gausspnt.RadiusArray = ['POINTS', field_name]
2395 if nb_components > 1:
2396 v_comp = get_vector_component(vector_mode)
2397 gausspnt.RadiusVectorComponent = v_comp
2398 gausspnt.RadiusTransferFunctionMode = 'Table'
2399 gausspnt.RadiusScalarRange = data_range
2400 gausspnt.RadiusUseScalarRange = 1
2401 if mult is not None:
2402 gausspnt.RadiusIsProportional = 1
2403 gausspnt.RadiusProportionalFactor = mult
2405 gausspnt.RadiusTransferFunctionEnabled = 0
2406 gausspnt.RadiusMode = 'Constant'
2407 gausspnt.RadiusArray = ['POINTS', 'Constant Radius']
2411 def StreamLinesOnField(proxy, entity, field_name, timestamp_nb,
2412 direction='BOTH', is_colored=False, color=None,
2413 vector_mode='Magnitude'):
2414 """Creates Stream Lines presentation on the given field.
2417 proxy: the pipeline object, containig data
2418 entity: the entity type from PrsTypeEnum
2419 field_name: the field name
2420 timestamp_nb: the number of time step (1, 2, ...)
2421 direction: the stream lines direction ('FORWARD', 'BACKWARD' or 'BOTH')
2422 is_colored: this option allows to color the presentation according to
2423 the corresponding data values. If False - the presentation will
2425 color: defines the presentation color as [R, G, B] triple. Taken into
2426 account only if is_colored is False.
2427 vector_mode: the mode of transformation of vector values
2428 into scalar values, applicable only if the field contains vector values.
2429 Possible modes: 'Magnitude', 'X', 'Y' or 'Z'.
2432 Stream Lines as representation object.
2435 proxy.UpdatePipeline()
2436 # We don't need mesh parts with no data on them
2437 if entity == EntityType.NODE:
2438 select_cells_with_data(proxy, on_points=[field_name])
2440 select_cells_with_data(proxy, on_cells=[field_name])
2443 nb_components = get_nb_components(proxy, entity, field_name)
2444 check_vector_mode(vector_mode, nb_components)
2447 time_value = get_time(proxy, timestamp_nb)
2450 pvs.GetRenderView().ViewTime = time_value
2451 pvs.UpdatePipeline(time_value, proxy)
2454 source = pvs.MergeBlocks(proxy)
2455 pvs.UpdatePipeline()
2457 # Cell data to point data
2458 if is_data_on_cells(proxy, field_name):
2459 cell_to_point = pvs.CellDatatoPointData(source)
2460 cell_to_point.PassCellData = 1
2461 pvs.UpdatePipeline()
2462 source = cell_to_point
2464 vector_array = field_name
2465 # If the given vector array has only 2 components, add the third one
2466 if nb_components == 2:
2467 calc = get_add_component_calc(source, EntityType.NODE, field_name)
2468 vector_array = calc.ResultArrayName
2469 pvs.UpdatePipeline()
2473 stream = pvs.StreamTracer(source)
2474 stream.SeedType = "Point Source"
2475 stream.Vectors = ['POINTS', vector_array]
2476 stream.IntegrationDirection = direction
2477 stream.IntegratorType = 'Runge-Kutta 2'
2478 stream.UpdatePipeline()
2480 # Get Stream Lines representation object
2481 if is_empty(stream):
2483 streamlines = pvs.GetRepresentation(stream)
2486 lookup_table = get_lookup_table(field_name, nb_components, vector_mode)
2488 # Set field range if necessary
2489 data_range = get_data_range(proxy, entity,
2490 field_name, vector_mode)
2491 if hasattr(lookup_table,"LockDataRange"):
2492 lookup_table.LockDataRange = 1
2493 elif hasattr(lookup_table,"LockScalarRange"):
2494 lookup_table.LockScalarRange = 1
2496 raise RuntimeError("Object %s has no 'LockDataRange' or 'LockScalarRange' attribute!"%(lookup_table))
2498 lookup_table.RGBPoints = [data_range[0], 0, 0, 1, data_range[1], 1, 0, 0]
2502 streamlines.ColorArrayName = (EntityType.get_pvtype(entity), field_name)
2504 streamlines.ColorArrayName = (None, '')
2506 streamlines.DiffuseColor = color
2508 streamlines.LookupTable = lookup_table
2511 add_scalar_bar(field_name, nb_components,
2512 vector_mode, lookup_table, time_value)
2517 def MeshOnEntity(proxy, mesh_name, entity):
2518 """Creates submesh of the entity type for the mesh.
2521 proxy -- the pipeline object, containig data
2522 mesh_name -- the full or short name of mesh field
2525 Submesh as representation object of the given source.
2528 proxy.UpdatePipeline()
2529 mesh_full_name = None
2530 aList = mesh_name.split('/')
2532 mesh_full_name = mesh_name
2534 mesh_full_name = find_mesh_full_name(proxy, mesh_name)
2535 if not mesh_full_name:
2536 raise RuntimeError, "The given mesh name was not found"
2537 # Select only the given mesh
2538 proxy.AllArrays = [mesh_full_name]
2539 proxy.UpdatePipeline()
2541 # Get representation object if the submesh is not empty
2543 if (proxy.GetDataInformation().GetNumberOfPoints() or
2544 proxy.GetDataInformation().GetNumberOfCells()):
2545 my_view = pvs.GetRenderView()
2546 prs = pvs.GetRepresentation(proxy, view=my_view)
2547 prs.ColorArrayName = (None, '')
2552 def MeshOnGroup(proxy, extrGroups, group_name):
2553 """Creates submesh on the group.
2556 proxy -- the pipeline object, containig data
2557 group_name -- the full group name
2558 extrGroups -- all extracted groups object
2561 Representation object of the given source with single group
2565 proxy.UpdatePipeline()
2566 # Deselect all groups
2567 extrGroups.AllGroups = []
2568 extrGroups.UpdatePipelineInformation()
2569 # Select only the group with the given name
2570 extrGroups.AllGroups = [group_name]
2571 extrGroups.UpdatePipelineInformation()
2573 # Get representation object if the submesh is not empty
2576 # Check if the group was set
2577 if len(extrGroups.AllGroups) == 1 and \
2578 extrGroups.AllGroups[0] == group_name:
2579 # Check if the submesh is not empty
2580 nb_points = proxy.GetDataInformation().GetNumberOfPoints()
2581 nb_cells = proxy.GetDataInformation().GetNumberOfCells()
2583 if nb_points or nb_cells:
2584 # prs = pvs.GetRepresentation(proxy)
2586 prs.ColorArrayName = ''
2592 def CreatePrsForFile(file_name, prs_types,
2593 picture_dir, picture_ext):
2594 """Build presentations of the given types for the file.
2596 Build presentations for all fields on all timestamps.
2599 file_name: full path to the MED file
2600 prs_types: the list of presentation types to build
2601 picture_dir: the directory path for saving snapshots
2602 picture_ext: graphics files extension (determines file type)
2606 print "Import " + file_name.split(os.sep)[-1] + "..."
2609 proxy = pvs.MEDReader(FileName=file_name)
2613 #proxy.UpdatePipeline()
2619 view = pvs.GetRenderView()
2620 time_value = get_time(proxy, 0)
2621 view.ViewTime = time_value
2622 pvs.UpdatePipeline(time=time_value, proxy=proxy)
2624 # Create required presentations for the proxy
2625 CreatePrsForProxy(proxy, view, prs_types,
2626 picture_dir, picture_ext)
2628 def CreatePrsForProxy(proxy, view, prs_types, picture_dir, picture_ext):
2629 """Build presentations of the given types for all fields of the proxy.
2631 Save snapshots in graphics files (type depends on the given extension).
2632 Stores the files in the given directory.
2635 proxy: the pipeline object, containig data
2636 view: the render view
2637 prs_types: the list of presentation types to build
2638 picture_dir: the directory path for saving snapshots
2639 picture_ext: graphics files extension (determines file type)
2642 proxy.UpdatePipeline()
2643 # List of the field names
2644 fields_info = proxy.GetProperty("FieldsTreeInfo")[::2]
2646 # Add path separator to the end of picture path if necessery
2647 if not picture_dir.endswith(os.sep):
2648 picture_dir += os.sep
2651 if PrsTypeEnum.MESH in prs_types:
2653 mesh_names = get_mesh_full_names(proxy)
2654 for mesh_name in mesh_names:
2655 # Build mesh field presentation
2656 print "Creating submesh for '" + get_field_short_name(mesh_name) + "' mesh... "
2657 prs = MeshOnEntity(proxy, mesh_name, None)
2663 # Construct image file name
2664 pic_name = picture_dir + get_field_short_name(mesh_name) + "." + picture_ext
2666 # Show and dump the presentation into a graphics file
2667 process_prs_for_test(prs, view, pic_name, False)
2669 # Create Mesh presentation. Build all groups.
2670 extGrp = pvs.ExtractGroup()
2671 extGrp.UpdatePipelineInformation()
2672 if if_possible(proxy, None, None, PrsTypeEnum.MESH, extGrp):
2673 for group in get_group_names(extGrp):
2674 print "Creating submesh on group " + get_group_short_name(group) + "... "
2675 prs = MeshOnGroup(proxy, extGrp, group)
2681 # Construct image file name
2682 pic_name = picture_dir + get_group_short_name(group) + "." + picture_ext
2684 # Show and dump the presentation into a graphics file
2685 process_prs_for_test(prs, view, pic_name, False)
2687 # Presentations on fields
2688 for field in fields_info:
2689 field_name = get_field_short_name(field)
2690 # Ignore mesh presentation
2691 if field_name == get_field_mesh_name(field):
2693 field_entity = get_field_entity(field)
2694 # Select only the current field:
2695 # necessary for getting the right timestamps
2696 proxy.AllArrays = [field]
2697 proxy.UpdatePipeline()
2700 timestamps = proxy.TimestepValues.GetData()
2702 for prs_type in prs_types:
2703 # Ignore mesh presentation
2704 if prs_type == PrsTypeEnum.MESH:
2707 # Get name of presentation type
2708 prs_name = PrsTypeEnum.get_name(prs_type)
2710 # Build the presentation if possible
2711 possible = if_possible(proxy, field_name,
2712 field_entity, prs_type)
2714 # Presentation type for graphics file name
2715 f_prs_type = prs_name.replace(' ', '').upper()
2717 for timestamp_nb in xrange(1, len(timestamps) + 1):
2718 time = timestamps[timestamp_nb - 1]
2720 scalar_range = get_data_range(proxy, field_entity,
2721 field_name, cut_off=True)
2722 # exclude time stamps with null lenght of scalar range
2723 if (scalar_range[0] == scalar_range[1]):
2725 print "Creating " + prs_name + " on " + field_name + ", time = " + str(time) + "... "
2727 prs = create_prs(prs_type, proxy,
2728 field_entity, field_name, timestamp_nb)
2730 """ This exception comes from get_nb_components(...) function.
2731 The reason of exception is an implementation of MEDReader
2732 activating the first leaf when reading MED file (refer to
2733 MEDFileFieldRepresentationTree::activateTheFirst() and
2734 MEDFileFieldRepresentationTree::getTheSingleActivated(...) methods).
2736 print "ValueError exception is catched"
2744 # Construct image file name
2745 pic_name = picture_dir + field_name + "_" + str(time) + "_" + f_prs_type + "." + picture_ext
2747 # Show and dump the presentation into a graphics file
2748 process_prs_for_test(prs, view, pic_name)