add script for coloured ply file conversion of the results

convert.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+#
+# Quick-and-dirty script which creates a coloured .ply file from a FE calculation.
+#   !! works with nodal values only !!
+#
+# the results can be viewed with "ctmviewer" (linux)
+#   (sudo apt install openctm-tools)
+#
+# usage:
+#   ./fossils.py models/dolicorhynchops/dolicorhynchops_10k.py
+#   ./convert.py
+#   ctmviewer philippe.ply
+
+import colorsys
+import vtk
+import gmsh
+import os.path
+import numpy as np
+
+import fossils
+fossils.setup_pythonpath()
+
+folder = os.path.join('workspace', \
+        'models_dolicorhynchops_dolicorhynchops_10k')
+meshfile = os.path.join(folder, 'mesh.msh')
+datafile = os.path.join(folder, 'smooth_stress_zz.msh') # <= values for colours
+
+# load mesh and result file into Gmsh
+gmsh.initialize()
+gmsh.merge(meshfile)
+gmsh.merge(datafile)
+
+# retrieve view data as numpy arrays
+views = gmsh.view.getTags()
+datatype, tags, data, \
+time, numComponents = gmsh.view.getHomogeneousModelData(views[0], 0)
+
+print(f'datatype = {datatype}')
+print(f'tags = {tags}')
+print(f'data = {data}')
+print(f'time = {time}')
+print(f'numComponents = {numComponents}')
+
+tmp_file = 'mesh.tmp.vtk'
+# write mesh to a file in VTK format
+gmsh.write(tmp_file)
+
+# read exported mesh
+reader = vtk.vtkUnstructuredGridReader()
+reader.SetFileName(tmp_file)
+reader.Update()
+ugrid = reader.GetOutput()
+
+# add a colour array to the mesh
+colors = vtk.vtkUnsignedCharArray()
+colors.SetNumberOfComponents(3)
+colors.SetName("Colors")
+ugrid.GetPointData().AddArray(colors)
+points = ugrid.GetPoints()
+npts = points.GetNumberOfPoints()
+# print(f'npts={npts}')
+# print(f'len(data)={len(data)}')
+vmin = np.min(data)
+vmax = np.max(data)
+# print(f'vmin={vmin}')
+# print(f'vmax={vmax}')
+rng = vmax-vmin
+
+for i, v in enumerate(data):
+    val = int((v-vmin)/rng*255) # compute "hue" value from 0 to 255
+    r, g, b = colorsys.hsv_to_rgb(val/255, 1.0, 1.0) # convert to "rainbow colours"
+    colors.InsertTuple3(i, r*255, g*255, b*255)
+
+# extract mesh boundary
+bfilter = vtk.vtkGeometryFilter()
+bfilter.SetInputConnection(reader.GetOutputPort())
+
+# export to .PLY
+writer = vtk.vtkPLYWriter()
+writer.SetFileName("philippe.ply")
+writer.SetArrayName("Colors")
+writer.SetInputConnection(bfilter.GetOutputPort())
+writer.Write()
+
+# clean tmp file
+os.remove(tmp_file)
+
+print('done.')