Maybe this snippets help new users. They are based on the scripts in the pytel folder.
The first thing we have to do is to import the matplotlib and the pytel modules. First we import some standard python modules
import sys
from os import path
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import tri
import matplotlib.cm as cm
Now we add the path to pytel to be able to import the modules. You will have to change the path (PathToPytel).
sys.path.append(path.join( path.dirname(sys.argv[0]), \
'/PathToPytel/opentelemac/v6p2/pytel'))
from utils.files import getFileContent
from parsers.parserSortie import getValueHistorySortie
from parsers.parserSELAFIN import SELAFIN,getValueHistorySLF,parseSLF,\
crossLocateMeshSLF,xyLocateMeshSLF,\
getValuePolylineSLF,subsetVariablesSLF
from parsers.parserStrings import parseArrayPaires
Now we should be able to read selafin files. I don't know all functions of the reader (slf object). Here we read the mesh, I assume that we have a triangular grid.
slf = SELAFIN(filename)
xmin = np.min(slf.MESHX)
xmax = np.max(slf.MESHX)
ymin = np.min(slf.MESHY)
ymax = np.max(slf.MESHY)
elements = np.dstack((slf.MESHX[slf.IKLE],slf.MESHY[slf.IKLE]))
mesh = np.array(slf.IKLE)
x = slf.MESHX
y = slf.MESHY
Now we can check the values and timesteps
for i,name in enumerate(slf.VARNAMES):
print "%i %s" %(i, name)
for i,time in enumerate(slf.tags["times"]):
print "%i %s" %(i, time)
At this point I was not sure. I think it is possible to get the values of a timestep i with getValues(i). Then the value can bee chosen by the number of the value (valueNumber).
values = slf.getVALUES(i)
ztri = values[valueNumber]
Now we are able to plot the value ztri and the mesh. First we make a new figure.
fig = plt.figure()
ax=fig.add_subplot(111)
There are different ways to plot the data. The best choice is to use the triplotcontourf function. We can define the levels of the plot.
zmintri = ztri.min()
zmaxtri = ztri.max()
nlevels = 15
levels = np.arange(zmintri,zmaxtri,(zmaxtri - zmintri)/nlevels)
triplotcontourf = plt.tricontourf(x,y,mesh,ztri,levels)
plt.colorbar()
Plotting the mesh can be done by defining polygons ( which is very slow) or with triplot (faster).
plt.triplot(x, y, mesh, 'k-',alpha= 0.2)
Velocities are plotted with quiver. Here we assume that the velocities are the first two entries of values. (See in the pytel scripts how the velocity can be normalized to improve the quality and get a nice plot)
u = values[0]
v = values[1]
cs = plt.quiver(x,y,u,v)
Now it is time to set the limits of the axis and plot our file.
ax.set_xlim(xmin,xmax)
ax.set_ylim(ymin,ymax)
ax.axis('equal')
plt.show()
plt.close()
Hope this helps new matplotlib users.
