Salome HOME
split zone tools - disable debtrace mess in release mode
[modules/hydro.git] / src / HYDROTools / interpolZ.py
index efe490a266b6685de2c808cff673992c3d7bd74c..93daf47ba89e217ad2e0f9da190e8466626f74db 100644 (file)
@@ -109,7 +109,7 @@ def createZfield1(fichierMaillage):
 
 # -----------------------------------------------------------------------------
 
-from MEDLoader import MEDLoader, MEDCouplingFieldDouble, ON_NODES, DataArrayDouble
+from MEDLoader import MEDLoader, MEDCouplingFieldDouble, ON_NODES, DataArrayDouble, MEDFileMesh
 
 def createZfield2(fichierMaillage):
   """
@@ -117,24 +117,24 @@ def createZfield2(fichierMaillage):
   Add a field on nodes, named "BOTTOM", of type double, containing z coordinates of nodes.
   createZfield2 is used after interpolZ. createZfield1 is base on MEDLoader interface.
   There is an alternate method based on Med file, equivalent (createZfield1).
-  The file <fichierMaillage>F.med produced by interpolz must exist, and is modified.
   fichierMaillage : 2D (x,y) med file produced by SMESH and used by interpolZ.
-  return <fichierMaillage>L.med : med file containing the field "BOTTOM"
+  return <fichierMaillage>F.med : med file containing the field "BOTTOM"
   """
   
   noms = string.split(fichierMaillage,'.')
   basename = string.join(noms[:-1], '.')
   fichierZMaillage = basename + 'Z.med'
-  fichierLMaillage = basename + 'L.med'
-  print fichierLMaillage
+  fichierFMaillage = basename + 'F.med'
+  print fichierFMaillage
 
   mymesh = MEDLoader.ReadUMeshFromFile(fichierZMaillage,0)
   fieldOnNodes=MEDCouplingFieldDouble.New(ON_NODES)
   fieldOnNodes.setName("BOTTOM")
   fieldOnNodes.setMesh(mymesh)
   fieldOnNodes.setArray(mymesh.getCoords()[:,2])
-
-  MEDLoader.WriteField(fichierLMaillage,fieldOnNodes,True)
+  mm=MEDFileMesh.New(fichierZMaillage)
+  mm.write(fichierFMaillage,2)
+  MEDLoader.WriteFieldUsingAlreadyWrittenMesh(fichierFMaillage,fieldOnNodes)
 
 # -----------------------------------------------------------------------------
 
@@ -156,7 +156,7 @@ from salome.smesh import smeshBuilder
 
 smesh = smeshBuilder.New(theStudy)
 
-def interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef):
+def interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef, interpolMethod = 0):
   """
   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.
@@ -166,6 +166,7 @@ def interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef):
   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
   """
@@ -190,6 +191,13 @@ def interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef):
   ([maillagePlat], status) = smesh.CreateMeshesFromMED(fichierMaillage)
   groups = maillagePlat.GetGroups()
   
+  grpns = [grp for grp in groups if grp.GetType() == SMESH.NODE]
+  if len(grpns) == 0:
+    print "Problem! There are no groups of nodes in the mesh!"
+    print "Please create at least the groups of nodes corresponding to each region of the HYDRO case" 
+    return {}
+    
+  
   for grp in groups:
     if grp.GetType() == SMESH.NODE:
       grpName = grp.GetName()
@@ -213,7 +221,7 @@ def interpolZ(nomCas, fichierMaillage, dicoGroupeRegion, zUndef):
       #print xyz
       vx.append(xyz[0])
       vy.append(xyz[1])
-    vz = cas.GetAltitudesForPoints( vx, vy, region )
+    vz = cas.GetAltitudesForPoints( vx, vy, region, interpolMethod )
     minz = min(vz)
     maxz = max(vz)
     statz[grp.GetName()] = (minz, maxz)