zUndef = 90
# --- Z interpolation on the bathymety/altimetry on the mesh nodes
-interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef)
-
-# --- add a field on nodes of type double with z values, named "BOTTOM"
-createZfield2(fichierMaillage)
+interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef, interpolMethod)
"""
-__revision__ = "V2.01"
+__revision__ = "V2.02"
# -----------------------------------------------------------------------------
# -----------------------------------------------------------------------------
# -----------------------------------------------------------------------------
-def interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef, interpolMethod=0, verbose=False):
+def interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef=90., interpolMethod=0, zUndefThreshold=-9000., xyzFile=False, verbose=False):
"""
interpolZ takes a 2D (x,y) mesh and calls the active instance of module HYDRO
to interpolate the bathymetry/altimetry on the mesh nodes, to produce the Z value of each node.
- interpolZ must be followed by createZfield1 or createZfield2.
- nomCas: Calculation Case Name in module HYDRO
- fichierMaillage: med file name produced by SMESH, corresponding to the HYDRO case
- dicoGroupeRegion: python dictionary giving the coorespondance of mesh groups to HYDRO regions.
- Key: face group name, value: region name in the HYDRO Case
- zUndef: Z value to use for nodes outside the regions (there must be none if the case is correct).
- interpolMethod: integer value, default 0 = nearest point on bathymetry, 1 = linear interpolation
- return <fichierMaillage>Z.med : med file with Z value on nodes
- return <fichierMaillage>F.med : an exact copy of <fichierMaillage>Z.med
+
+ In:
+ nomCas: Calculation Case Name in module HYDRO
+ fichierMaillage: med file name produced by SMESH, corresponding to the HYDRO case
+ dicoGroupeRegion: python dictionary giving the correspondance of mesh groups to HYDRO regions.
+ Key: face group name, value: region name in the HYDRO Case
+ zUndef: Z value to use for nodes outside the regions (there must be none if the case is correct).
+ default value is 90.
+ interpolMethod: integer value
+ 0 = nearest point on bathymetry (default)
+ 1 = linear interpolation
+ zUndef: Z value to use for nodes outside the regions (there must be none if the case is correct).
+ zUndefThreshold: if the interpolated z is under this value, it is replaced by zUndef
+ default value is -9000.
+ xyzFile: True/False to write an ascii file with xyz for every node. Default is False.
+ Out:
+ statz: statistique for z
+ Key: face group name, value: (minz, maxz, meanz, stdz, v05z, v95z)
+ Out:
+ return <fichierMaillage>F.med : med file with Z value on nodes and in a field "BOTTOM"
"""
statz = dict()
erreur = 0
message = ""
while not erreur:
+
+ if verbose:
+ print "interpolMethod: %d" % interpolMethod
+ print "zUndef:", zUndef
+ print "zUndefThreshold: %f" % zUndefThreshold
doc = HYDROPy.HYDROData_Document.Document(theStudyId)
cas = doc.FindObjectByName(nomCas)
custom_inter = MyInterpolator()
basename = fichierMaillage[:-4]
- fichierZMaillage = basename + 'Z.med'
fichierFMaillage = basename + 'F.med'
- fichierFonds = basename + '.xyz'
- print "fichierMaillage :", fichierMaillage
- print "fichierZMaillage :", fichierZMaillage
- print "fichierFMaillage :", fichierFMaillage
- print "fichierFonds :", fichierFonds
print "dicoGroupeRegion =", dicoGroupeRegion
+ print "fichierMaillage =", fichierMaillage
+ print "fichierFMaillage =", fichierFMaillage
+ if xyzFile:
+ fichierFonds = basename + '.xyz'
+ print "fichierFonds =", fichierFonds
#
# 1. Reads the mesh and gets the number of nodes
#
message += "This group does"
else:
message += "That %d groups do" % nb_pb
- message += " not belongs to the mesh\n"
- message += "Please check the names of the groups of faces corresponding to each region of the HYDRO case"
+ message += " not belongs to the mesh.\n"
+ message += "Please check the names of the group(s) of faces corresponding to each region of the HYDRO case"
erreur = 2
break
#
-# 3. Get the information about the nodes
+# 3. Gets the information about the nodes
#
nbnodes = meshMEDFileRead.getNumberOfNodes()
if verbose:
#print coords
#print coords[0,0]
#print coords[0,1]
+#
+# 4. Exploration of every group of faces
#
tb_aux = np.zeros(nbnodes, dtype=np.bool)
- bathy = np.zeros(nbnodes, dtype=np.float)
#
-# 4. Exploration of every group of faces
+ bathy = np.zeros(nbnodes, dtype=np.float)
+ bathy.fill(zUndef)
#
for gr_face_name in l_gr_faces:
#
-# 4.1. Region
+# 4.1. Region connected to the group
#
nomreg = dicoGroupeRegion[gr_face_name]
line = "------- Region: '" + nomreg + "'"
print "\t. Number of cells: %d" % nbr_cells
#
# 4.3. Nodes of the meshes of the group
+# Every node is flagged in tb_aux
#
tb_aux.fill(False)
for id_elem in range(nbr_cells):
np_aux = tb_aux.nonzero()
if len(np_aux[0]):
if verbose:
- print "\t. Number of nodes: %d" % len(np_aux[0])
+ print "\t. Number of nodes for this group: %d" % len(np_aux[0])
#print "np_aux:", np_aux
#
# 4.4. Interpolation over the nodes of the meshes of the group
ligne += ", v95z: %f" % v95z
print ligne
#
-# 4.5. Storage
+# 4.5. Storage of the z and of the statistics for this region
#
statz[gr_face_name] = (minz, maxz, meanz, stdz, v05z, v95z)
+#
for iaux, nodeId in enumerate(np_aux[0]):
bathy[nodeId] = vz[iaux]
#
-# 5. Minimum
+# 5. Minimum:
+# If the value is lower than a threshold, an "undefined" valeur is set
#
#print "bathy :\n", bathy
- np_aux_z = (bathy < -9000.).nonzero()
+ np_aux_z = (bathy < zUndefThreshold).nonzero()
if verbose:
print ".. Number of nodes below the minimum: %d" % len(np_aux_z[0])
if len(np_aux_z[0]):
#
# 6. xyz file
#
- if verbose:
- print ".. Ecriture du champ de bathymétrie sur le fichier :\n", fichierFonds
- fo = open(fichierFonds, 'w')
- for nodeId in range(nbnodes):
- #maillagePlat.MoveNode(nodeId, x, y, z)
- line = "%10.2f %10.2f %10.2f\n" % (coords[nodeId, 0], coords[nodeId, 1], bathy[nodeId])
- fo.write(line)
- fo.close()
+ if xyzFile:
#
-# 7 Modification of the z coordinates
+ if verbose:
+ print ".. Ecriture du champ de bathymétrie sur le fichier :\n", fichierFonds
+ fo = open(fichierFonds, 'w')
+ for nodeId in range(nbnodes):
+ line = "%10.2f %10.2f %10.2f\n" % (coords[nodeId, 0], coords[nodeId, 1], bathy[nodeId])
+ fo.write(line)
+ fo.close()
+#
+# 7. Final MED file
+# 7.1. Modification of the z coordinates
#
bathy_dd = mc.DataArrayDouble(np.asfarray(bathy, dtype='float'))
coords3D = DataArrayDouble.Meld([coords, bathy_dd])
#print "coords3D\n", coords3D
meshMEDFileRead.setCoords(coords3D)
+#
+# 7.2. Writes the 3D mesh
#
if verbose:
- print ".. Ecriture du maillage 3D sur le fichier :\n", fichierZMaillage
- meshMEDFileRead.write(fichierZMaillage, 2)
-
- if verbose:
- print ".. Ecriture du maillage 3D avec le champ BOTTOM sur le fichier :\n", fichierFMaillage
+ print ".. Ecriture du maillage 3D sur le fichier :\n", fichierFMaillage
meshMEDFileRead.write(fichierFMaillage, 2)
-
+#
+# 7.3. Writes the field
+#
+ if verbose:
+ print ".. Ecriture du champ BOTTOM"
+ fieldOnNodes=ml.MEDCouplingFieldDouble(ml.ON_NODES, ml.ONE_TIME)
+ fieldOnNodes.setName("BOTTOM")
+ fieldOnNodes.setMesh(meshMEDFileRead.getMeshAtLevel(0))
+ fieldOnNodes.setArray(bathy_dd)
+#
+ fMEDFile_ch_d = ml.MEDFileField1TS()
+ fMEDFile_ch_d.setFieldNoProfileSBT(fieldOnNodes)
+ fMEDFile_ch_d.write(fichierFMaillage, 0)
+#
break
if erreur:
