src/HYDROTools/interpolS.py

   1 # -*- coding: utf-8 -*-
   2
   3 import os
   4 import sys
   5 import shutil
   6
   7 import salome
   8
   9 from med import medfile
  10 from med import medmesh
  11 from med import medfield
  12 from med import medenum
  13
  14 import HYDROPy
  15
  16 # Get current study id
  17 salome.salome_init()
  18 theStudy = salome.myStudy
  19 #theStudyId = salome.myStudy._get_StudyId()
  20
  21 """Bad parameters exception"""
  22 class BadParamsError(ValueError):
  23   pass
  24
  25 """
  26 Creates MED file with scalar field associated with the mesh of the calculation case.
  27 The scalar value associated with the mesh node is the Strickler coefficient of the land cover containing the node.
  28 case_name: calculation case name in the study
  29 med_file_name: path to input MED file with mesh on nodes
  30 output_file_name: path to output MED file with
  31 med_field_name: field name
  32 """
  33 def assignStrickler(case_name, med_file_name, output_file_name, med_field_name='FRICTION'):
  34   # Check calculation case
  35   doc = HYDROPy.HYDROData_Document.Document()
  36   case = doc.FindObjectByName(case_name)
  37   if case is None:
  38     raise BadParamsError("Calculation case '%s' not found" % case_name)
  39
  40   # Check input MED file
  41   if not os.path.exists(med_file_name):
  42     raise BadParamsError("Input file '%s' not exists" % med_file_name)
  43
  44   # Copy input file to the output path, if the output path is empty - use the input path for the output
  45   file_path = med_file_name
  46   if output_file_name and (output_file_name != file_path):
  47     shutil.copyfile(med_file_name, output_file_name)
  48     file_path = output_file_name
  49
  50   # Open input MED file
  51   fid = medfile.MEDfileOpen(file_path, medenum.MED_ACC_RDEXT)
  52
  53   # Get mesh info
  54   mesh_name, sdim, mdim, meshtype, desc, dtunit, sort, nstep,  repere, axisname, axisunit = medmesh.MEDmeshInfo(fid, 1)
  55
  56   nb_nodes, chgt, trsf = medmesh.MEDmeshnEntity(fid, mesh_name, medenum.MED_NO_DT, medenum.MED_NO_IT,
  57                                                 medenum.MED_NODE, medenum.MED_NONE,
  58                                                 medenum.MED_COORDINATE, medenum.MED_NO_CMODE)
  59
  60
  61   # Create field
  62   if nb_nodes > 0:
  63     # Get node coordinates
  64     coords = medfile.MEDFLOAT(nb_nodes*sdim)
  65     medmesh.MEDmeshNodeCoordinateRd(fid, mesh_name, medenum.MED_NO_DT, medenum.MED_NO_IT,
  66                                     medenum.MED_FULL_INTERLACE, coords)
  67
  68     # Get list of Strickler coefficient values for the nodes
  69     values = []
  70     x_coords = coords[0::sdim]
  71     y_coords = coords[1::sdim]
  72     values = case.GetStricklerCoefficientForPoints(x_coords, y_coords, 0.0, True)
  73
  74     # Write the values to the field
  75     if len(values) > 0:
  76       values = medfile.MEDFLOAT(values)
  77
  78       comp  = "strickler coeff "
  79       medfield.MEDfieldCr(fid, med_field_name, medfile.MED_FLOAT64,
  80                           1, comp, '', '', mesh_name)
  81
  82       medfield.MEDfieldValueWr(fid, med_field_name, medenum.MED_NO_DT, medenum.MED_NO_IT, 0.0,
  83                                medenum.MED_NODE, medenum.MED_NONE,
  84                                medenum.MED_FULL_INTERLACE, medenum.MED_ALL_CONSTITUENT, len(values), values)
  85
  86       print("MED field '%s' on %s nodes has been created." % (med_field_name, nb_nodes))
  87
  88   # close MED file
  89   medfile.MEDfileClose(fid)