src/MEDMEM_SWIG/medutilities.py

   1 # -*- coding: iso-8859-1 -*-
   2 # --
   3 # Copyright (C) 2009-2013  CEA/DEN, EDF R&D, OPEN CASCADE
   4 #
   5 # This library is free software; you can redistribute it and/or
   6 # modify it under the terms of the GNU Lesser General Public
   7 # License as published by the Free Software Foundation; either
   8 # version 2.1 of the License.
   9 #
  10 # This library is distributed in the hope that it will be useful,
  11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 # Lesser General Public License for more details.
  14 #
  15 # You should have received a copy of the GNU Lesser General Public
  16 # License along with this library; if not, write to the Free Software
  17 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  18 #
  19 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  20 #
  21 # Author : Erwan ADAM (CEA)
  22 # --
  23
  24 from libMEDMEM_Swig import *
  25
  26 def my_remove(f):
  27     from os import remove
  28     try:
  29         remove(f)
  30     except OSError:
  31         pass
  32     return
  33
  34 def sauv2med(*argv):
  35     argv = list(argv)
  36     # argv = argv[1:]
  37     for arg in argv:
  38         convert(arg, "GIBI", "MED")
  39         pass
  40     return
  41
  42 def med2sauv(*argv):
  43     argv = list(argv)
  44     # argv = argv[1:]
  45     format = 1
  46     for arg in argv[:]:
  47         if arg.find('--format') == 0:
  48             argv.remove(arg)
  49             try:
  50                 value = arg.split("=")[1]
  51             except IndexError:
  52                 usage(1)
  53                 pass
  54             try:
  55                 value = int(value)
  56             except ValueError:
  57                 usage(1)
  58                 pass
  59             format = value
  60             pass
  61         pass
  62     for arg in argv:
  63         convert(arg, "MED", "GIBI", format)
  64         pass
  65     return
  66
  67 def convert(file_in, driver_in, driver_out, format=1, file_out=None):
  68     #
  69     print file_in
  70     #
  71     if file_out is None:
  72         file_out = file_in
  73 ##        if file_out.find('.') != -1:
  74 ##            suffix = file_in.split('.')[-1]
  75 ##            if driver_in == "GIBI":
  76 ##                test = "sauv"
  77 ##            else:
  78 ##                test = "med"
  79 ##                pass
  80 ##            if len(suffix) >= len(test):
  81 ##                suffix = suffix[:len(test)]
  82 ##                suffix = suffix.lower()
  83 ##                if suffix == test:
  84 ##                    file_out = '.'.join(file_in.split('.')[:-1])
  85 ##                    pass
  86 ##                pass
  87 ##            pass
  88         if driver_out == "GIBI":
  89             file_out += ".sauv"
  90         elif driver_out == "MED":
  91             file_out += ".med"
  92         else:
  93             msg = "Driver out %s is unknown"%(driver_out)
  94             raise NotImplementedError(msg)
  95         pass
  96     print file_out
  97     #
  98     meshes = []
  99     fields = []
 100     if driver_in == "GIBI":
 101         driver = GIBI_MED_RDONLY_DRIVER(file_in)
 102         fields = driver.read()
 103         mesh = driver.getMesh()
 104         if mesh:
 105             meshes.append( mesh )
 106     elif driver_in == "MED":
 107         med = MEDFILEBROWSER(file_in)
 108         for mesh_name in med.getMeshNames():
 109             if med.isStructuredMesh( mesh_name ):
 110                 mesh = GRID( MED_DRIVER, file_in, mesh_name )
 111             else:
 112                 mesh = MESH( MED_DRIVER, file_in, mesh_name )
 113             meshes.append( mesh )
 114         for field_name in med.getFieldNames():
 115             mesh_name = med.getMeshName( field_name )
 116             mesh = 0
 117             for m in meshes:
 118                 if m.getName() == mesh_name:
 119                     mesh = m; break
 120             for dtit in med.getFieldIteration( field_name ):
 121                 if med.getFieldType( field_name ) == MED_REEL64:
 122                     field = FIELDDOUBLE(MED_DRIVER, file_in, field_name, dtit.dt, dtit.it, mesh )
 123                 else:
 124                     field = FIELDINT(MED_DRIVER, file_in, field_name, dtit.dt, dtit.it, mesh )
 125                 fields.append( field )
 126     else:
 127         msg = "Driver in %s is unknown"%(driver_in)
 128         raise NotImplementedError(msg)
 129     #
 130     my_remove(file_out)
 131     #
 132     if driver_out == "GIBI":
 133         mesh = meshes[0]
 134         mesh_dim = mesh.getSpaceDimension()
 135         if format == 0:
 136             file_out = file_out+'__format__'
 137             my_remove(file_out)
 138             pass
 139         if fields:
 140             driver = GIBI_MED_WRONLY_DRIVER(file_out, fields, mesh)
 141         else:
 142             driver = GIBI_MESH_WRONLY_DRIVER(file_out, mesh)
 143         driver.open()
 144         driver.write()
 145         driver.close()
 146         #
 147         if mesh_dim >= 3:
 148             from sys import platform
 149             if platform in ["win32"]:
 150                 f = open(file_out)
 151                 content = f.read()
 152                 f.close()
 153                 content = content.replace("IFOUR  -1", "IFOUR   2")
 154                 content = content.replace("IFOMOD  -1", "IFOMOD   2")
 155                 f = open(file_out, "w")
 156                 f.write(content)
 157                 f.close()
 158             else:
 159                 cmd  = "sed"
 160                 cmd += ' -e "s/IFOUR  -1/IFOUR   2/g"'
 161                 cmd += ' -e "s/IFOMOD  -1/IFOMOD   2/g"'
 162                 # cmd += ' -e "s/IECHO   1/IECHO   0/g"'
 163                 cmd += ' %s > .dummy'%(file_out)
 164                 cmd += ' && '
 165                 cmd += ' mv -f .dummy %s'%(file_out)
 166                 from os import system
 167                 system(cmd)
 168                 pass
 169             pass
 170         #
 171         if format == 0:
 172             from castemlauncher import CastemLauncher
 173             dgibi_stream  = "\n"
 174             dgibi_stream += "OPTI REST FORMAT '%s' ;\n"%(file_out)
 175             dgibi_stream += "REST FORMAT;\n"
 176             file_out = file_out.replace('__format__', '')
 177             dgibi_stream += "OPTI SAUV '%s' ;\n"%(file_out)
 178             dgibi_stream += "SAUV ;\n"
 179             cl = CastemLauncher(dgibi_stream)
 180             cl.addTmpFiles(file_out+'__format__', "UTILNOTI", "UTILPROC")
 181             cl.run()
 182             pass
 183         return
 184     #
 185     for mesh in meshes:
 186         mesh.write(MED_DRIVER, file_out)
 187     for field in fields:
 188         typedField = field.castToTypedField();
 189         typedField.write(MED_DRIVER, file_out)
 190     #
 191     return
 192
 193 def avs2med_one_file(file_in, file_out, mesh_name, field_name):
 194     """
 195     Convert an ucd avs inp file into a med file
 196     inp Specifications can be found at :
 197     http://people.scs.fsu.edu/~burkardt/data/ucd/ucd.html
 198     http://help.avs.com/Express/doc/help/reference/dvmac/UCD_Form.htm
 199
 200     """
 201     my_remove(file_out)
 202     #
 203     meshing = MESHING()
 204     meshing.setName(mesh_name)
 205     #
 206     f = open(file_in)
 207     lines = f.readlines()
 208     f.close()
 209     nb_lines = len(lines)
 210     # ----
 211     # Skip the comments
 212     # ----
 213     while 1:
 214         l = lines[0]
 215         if l[0] != "#":
 216             break
 217         lines = lines[1:]
 218         pass
 219     # ----
 220     headers = [ int(i) for i in lines[0].split() ]
 221     lines = lines[1:]
 222     number_of_nodes = headers[0]
 223     number_of_cells = headers[1]
 224     number_of_nodes_data = headers[2]
 225     number_of_cells_data = headers[3]
 226     number_of_whole_data = headers[3]
 227     # ----------
 228     # Nodes
 229     # ----------
 230     nodes = lines[:number_of_nodes]
 231     lines = lines[number_of_nodes:]
 232     nodes = [ " ".join(l.split()[1:]) for l in nodes ]
 233     nodes = " ".join(nodes)
 234     nodes = [ float(v) for v in nodes.split() ]
 235     # --------
 236     # Space dimension
 237     # --------
 238     nodes_min = min(nodes)
 239     nodes_max = max(nodes)
 240     epsilon = 1.0e-5 * (nodes_max - nodes_min)
 241     nodes_z = [ nodes[i] for i in range(2, 3 * number_of_nodes, 3) ]
 242     nodes_z_min = min(nodes_z)
 243     nodes_z_max = max(nodes_z)
 244     space_dimension = 2
 245     if abs(nodes_z_max) > epsilon:
 246         space_dimension = 3
 247         pass
 248     if abs(nodes_z_min) > epsilon:
 249         space_dimension = 3
 250         pass
 251     #
 252     if space_dimension == 2:
 253         l = range(0, 3 * number_of_nodes, 3) + range(1, 3 * number_of_nodes, 3)
 254         l.sort()
 255         nodes = [ nodes[i] for i in l ]
 256         pass
 257     #
 258     meshing.setCoordinates(space_dimension, number_of_nodes, nodes, "CARTESIAN", MED_FULL_INTERLACE)
 259     # ----------
 260     # Cells
 261     # ----------
 262     cells = lines[:number_of_cells]
 263     lines = lines[number_of_cells:]
 264     cells = [ c.split() for c in cells ]
 265     #
 266     type2connectivity = {}
 267     for c in cells:
 268         cell_id = int(c[0])
 269         cell_mat = int(c[1])
 270         cell_type = c[2]
 271         cell_connectivity = c[3:]
 272         #
 273         avs_type = cell_type
 274         if 0:
 275             pass
 276         elif avs_type == 'tri':
 277             entity = MED_CELL
 278             nb2medtype = {
 279                 3: MED_TRIA3,
 280                 }
 281         elif avs_type == 'quad':
 282             entity = MED_CELL
 283             nb2medtype = {
 284                 4: MED_QUAD4,
 285                 8: MED_QUAD8,
 286                 }
 287         elif avs_type == 'tet':
 288             entity = MED_CELL
 289             nb2medtype = {
 290                 4: MED_TETRA4,
 291                 }
 292         elif avs_type == 'hex':
 293             entity = MED_CELL
 294             nb2medtype = {
 295                 8: MED_HEXA8,
 296                 }
 297         elif avs_type == 'pyr':
 298             entity = MED_CELL
 299             nb2medtype = {
 300                 5: MED_PYRA5,
 301                 }
 302         elif avs_type == 'prism':
 303             entity = MED_CELL
 304             nb2medtype = {
 305                 6: MED_PENTA6,
 306                 }
 307         else:
 308             raise Exception("unknown avs_type : %s"%(avs_type))
 309         #
 310         med_type = nb2medtype[len(cell_connectivity)]
 311         #
 312         try:
 313             d = type2connectivity[entity]
 314         except:
 315             type2connectivity[entity] = {}
 316             d = type2connectivity[entity]
 317             pass
 318         #
 319         try:
 320             l = d[med_type]
 321         except:
 322             d[med_type] = []
 323             l = d[med_type]
 324             pass
 325         #
 326         cell_connectivity = ' '.join(cell_connectivity)
 327         l.append(cell_connectivity)
 328         #
 329         pass
 330     #
 331     mesh_dimension = space_dimension
 332     #
 333     for entity, d in type2connectivity.items():
 334         meshing.setNumberOfTypes(len(d.keys()), entity)
 335         meshing.setTypes(d.keys(), entity)
 336         meshing.setNumberOfElements([len(v) for v in d.values()], entity)
 337         for key, value in d.items():
 338             connectivity = value
 339             connectivity = " ".join(connectivity)
 340             connectivity = [ int(v) for v in connectivity.split() ]
 341             meshing.setConnectivity(connectivity, entity, key)
 342             pass
 343         pass
 344     # -----------
 345     meshing.write(meshing.addDriver(MED_DRIVER,file_out,meshing.getName()))
 346     # -----------
 347     flag = -1
 348     for n in [number_of_nodes_data, number_of_cells_data]:
 349         flag += 1
 350         if n == 0:
 351             continue
 352         # -----------
 353         header = [ int(v) for v in lines[0].split() ]
 354         lines = lines[1:]
 355         nb_of_components = header[0]
 356         nb_of_values_by_component = header[1:]
 357         if len(nb_of_values_by_component) != nb_of_components:
 358             msg = "Error at line %d"%(nb_lines - len(lines))
 359             raise Exception(msg)
 360         #
 361         titles_by_component = []
 362         for i in range(nb_of_components):
 363             l = lines[0]
 364             lines = lines[1:]
 365             label = l.split(',')[0]
 366             unit = l[len(label)+1:]   # The +1 is for the comma
 367             label = label.strip()
 368             unit = unit.strip()
 369             labels_by_value = nb_of_values_by_component[i]*['']
 370             l = [label, unit, labels_by_value]
 371             titles_by_component.append(l)
 372             pass
 373         if nb_of_components > 1:
 374             if nb_of_values_by_component == nb_of_components * [1]:
 375                 nb_of_values_by_component = [ nb_of_components ]
 376                 nb_of_components = 1
 377                 if flag == 0:
 378                     if number_of_cells_data:
 379                         name = field_name + "_on_nodes"
 380                         pass
 381                     pass
 382                 else:
 383                     if number_of_nodes_data:
 384                         name = field_name + "_on_cells"
 385                         pass
 386                     pass
 387                 titles_by_component = [ [name, "", [t[0] for t in titles_by_component]] ]
 388                 pass
 389             pass
 390         if flag == 0:
 391             nb = number_of_nodes
 392         else:
 393             nb = number_of_cells
 394             pass
 395         values = lines[:nb]
 396         lines = lines[nb:]
 397         vals = []
 398         imin = 1
 399         for i in range(nb_of_components):
 400             imax = imin + nb_of_values_by_component[i]
 401             vals.append([ l.split()[imin:imax] for l in values ])
 402             imin = imax
 403             pass
 404         values = vals
 405         # ----------
 406         if flag == 0:
 407             support = SUPPORT(meshing, "support_on_nodes", MED_NODE)
 408         else:
 409             support = SUPPORT(meshing, "support_on_cells", MED_CELL)
 410             pass
 411         for i in range(nb_of_components):
 412             nb = nb_of_values_by_component[i]
 413             field = FIELDDOUBLE(support, nb)
 414             field.setName(titles_by_component[i][0])
 415             # field.setIterationNumber(-1)
 416             # field.setOrderNumber(-1)
 417             # field.setTime(0.0)
 418             for n in range(nb):
 419                 name = titles_by_component[i][2][n]
 420                 if name:
 421                     field.setComponentName(n+1,name)
 422                     pass
 423                 pass
 424             # ---------------
 425             value = values[i]
 426             value = [ " ".join(l) for l in value ]
 427             value = " ".join(value)
 428             value = [ float(v) for v in value.split() ]
 429             field.setValue(value)
 430             # -----
 431             field.write(field.addDriver(MED_DRIVER,file_out,field.getName()))
 432             # -------
 433             pass
 434         #
 435         pass
 436     return
 437
 438 def avs2med(*argv):
 439     argv = list(argv)
 440     #
 441     mesh_name = "mesh"
 442     field_name = "field"
 443     #
 444     for arg in argv[:]:
 445         if ( arg.find('--mesh_name') == 0 ) or ( arg.find('--mesh-name') == 0 ) :
 446             argv.remove(arg)
 447             try:
 448                 mesh_name = arg.split("=")[1]
 449             except IndexError:
 450                 usage(1)
 451                 pass
 452             pass
 453         if ( arg.find('--field_name') == 0 ) or ( arg.find('--field-name') == 0 ) :
 454             argv.remove(arg)
 455             try:
 456                 field_name = arg.split("=")[1]
 457             except IndexError:
 458                 usage(1)
 459                 pass
 460             pass
 461         pass
 462     #
 463     for arg in argv:
 464         avs2med_one_file(arg, arg + ".med", mesh_name, field_name)
 465         pass
 466     #
 467     return
 468