+++ /dev/null
-# -*- coding: iso-8859-1 -*-\r
-# Copyright (C) 2007-2016 CEA/DEN, EDF R&D\r
-#\r
-# This library is free software; you can redistribute it and/or\r
-# modify it under the terms of the GNU Lesser General Public\r
-# License as published by the Free Software Foundation; either\r
-# version 2.1 of the License, or (at your option) any later version.\r
-#\r
-# This library is distributed in the hope that it will be useful,\r
-# but WITHOUT ANY WARRANTY; without even the implied warranty of\r
-# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU\r
-# Lesser General Public License for more details.\r
-#\r
-# You should have received a copy of the GNU Lesser General Public\r
-# License along with this library; if not, write to the Free Software\r
-# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA\r
-#\r
-# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com\r
-#\r
-# Author : Anthony GEAY (CEA/DEN/DM2S/STMF)\r
-\r
-from MEDLoader import *\r
-\r
-class PVDReader:\r
- @classmethod\r
- def New(cls,fileName):\r
- """ Static constructor. """\r
- return PVDReader(fileName)\r
- pass\r
-\r
- def __init__(self,fileName):\r
- self._fileName=fileName\r
- pass\r
-\r
- def loadTopInfo(self):\r
- fd=open(self._fileName,"r")\r
- return self.__parseXML(fd)\r
-\r
- def __parseXML(self,fd):\r
- import xml.sax\r
- class PVD_SAX_Reader(xml.sax.ContentHandler):\r
- def __init__(self):\r
- self._tsteps=[]\r
- pass\r
- def startElement(self,name,attrs):\r
- if name=="VTKFile":\r
- if attrs["type"]!="Collection":\r
- raise Exception("Mismatch between reader (PVD) type and file content !")\r
- return\r
- if name=="DataSet":\r
- self._tsteps.append((float(attrs["timestep"]),str(attrs["file"])))\r
- return\r
- pass\r
- pass\r
- rd=PVD_SAX_Reader()\r
- parser=xml.sax.make_parser()\r
- parser.setContentHandler(rd)\r
- parser.parse(fd)\r
- return rd\r
- pass\r
-\r
-class PVTUReader:\r
- @classmethod\r
- def New(cls,fileName):\r
- """ Static constructor. """\r
- return PVTUReader(fileName)\r
- pass\r
-\r
- def __init__(self,fileName):\r
- self._fileName=fileName\r
- pass\r
-\r
- def loadParaInfo(self):\r
- fd=open(self._fileName,"r")\r
- return self.__parseXML(fd)\r
-\r
- def __parseXML(self,fd):\r
- import xml.sax\r
- class PVTU_SAX_Reader(xml.sax.ContentHandler):\r
- def __init__(self):\r
- self._data_array={2:self.DAPointData,3:self.DACellData}\r
- self._node_fields=[]\r
- self._cell_fields=[]\r
- self._pfiles=[]\r
- self._tmp=None\r
- pass\r
- def DAPointData(self,attrs):\r
- self._node_fields.append((str(attrs["Name"]),int(attrs["NumberOfComponents"])))\r
- pass\r
- def DACellData(self,attrs):\r
- self._cell_fields.append((str(attrs["Name"]),int(attrs["NumberOfComponents"])))\r
- pass\r
- def startElement(self,name,attrs):\r
- if name=="VTKFile":\r
- if attrs["type"]!="PUnstructuredGrid":\r
- raise Exception("Mismatch between reader (PVTU) type and file content !")\r
- return\r
- if name=="Piece":\r
- self._pfiles.append(str(attrs["Source"]))\r
- return\r
- if name=="PPointData":\r
- self._tmp=2\r
- return\r
- if name=="PCellData":\r
- self._tmp=3\r
- return\r
- if name=="PDataArray":\r
- if self._tmp in self._data_array:\r
- self._data_array[self._tmp](attrs)\r
- pass\r
- return\r
- pass\r
- pass\r
- rd=PVTU_SAX_Reader()\r
- parser=xml.sax.make_parser()\r
- parser.setContentHandler(rd)\r
- parser.parse(fd)\r
- return rd\r
- pass\r
-\r
-class VTURawReader:\r
- """ Converting a VTU file in raw mode into the MED format.\r
- Warning: VTU file must be write in "Appended" mode, and in "Binary" format.\r
- """\r
- VTKTypes_2_MC=[-1,0,-1,1,33,NORM_TRI3 ,-1,5,-1,4,14,-1,NORM_HEXA8,16,15,-1,22,-1,-1,-1,-1,2,6,8,20,30,25,23,9,27,-1,-1,-1,-1,7,-1,-1,-1,-1,-1,-1,-1,31]\r
-\r
- class NormalException(Exception):\r
- def __init__(self,lineNb):\r
- Exception.__init__(self)\r
- self._line_nb=lineNb\r
- def getLineNb(self):\r
- return self._line_nb\r
- pass\r
-\r
- class NotRawVTUException(Exception):\r
- pass\r
-\r
- def loadInMEDFileDS(self):\r
- import numpy as np\r
- fd=open(self._fileName,"r")\r
- ref,rd=self.__parseXML(fd)\r
- #\r
- ret=MEDFileData()\r
- ms=MEDFileMeshes() ; ret.setMeshes(ms)\r
- fs=MEDFileFields() ; ret.setFields(fs)\r
- #\r
- types=np.memmap(fd,dtype=rd._type_types,mode='r',offset=ref+rd._off_types,shape=(rd._nb_cells,))\r
- types=self.__swapIfNecessary(rd._bo,types)\r
- # mesh dimension detection\r
- types2=types.copy() ; types2.sort() ; \r
- types2 =np.unique(types2)\r
- # Get first valid type\r
- meshDim = -1\r
- for i, typ in enumerate(types2):\r
- if self.VTKTypes_2_MC[typ] != -1:\r
- meshDim=MEDCouplingMesh.GetDimensionOfGeometricType(self.VTKTypes_2_MC[typ])\r
- if meshDim == -1:\r
- raise Exception("Could not find a valid cell type in the mesh !")\r
- if i > 0:\r
- print("WARNING: some invalid/incompatible cell types were detected - trying to have them as polygons ...")\r
- for typ in types2:\r
- if self.VTKTypes_2_MC[typ] != -1:\r
- md=MEDCouplingMesh.GetDimensionOfGeometricType(self.VTKTypes_2_MC[typ])\r
- if md!=meshDim:\r
- raise Exception("MultiLevel umeshes not managed yet !")\r
- else:\r
- print("WARNING: invalid/incompatible cell type detected: VTK type number: %d" % typ)\r
- m=MEDCouplingUMesh("mesh",meshDim)\r
- # coordinates\r
- coo=np.memmap(fd,dtype=rd._type_coords,mode='r',offset=ref+rd._off_coords,shape=(rd._nb_nodes*rd._space_dim,))\r
- coo=self.__swapIfNecessary(rd._bo,coo) ; coo=DataArrayDouble(np.array(coo,dtype='float64')) ; coo.rearrange(rd._space_dim)\r
- m.setCoords(coo)\r
- # connectivity\r
- offsets=np.memmap(fd,dtype=rd._type_off,mode='r',offset=ref+rd._off_off,shape=(rd._nb_cells,))\r
- offsets=self.__swapIfNecessary(rd._bo,offsets) ; connLgth=offsets[-1] ; offsets2=DataArrayIdType(rd._nb_cells+1) ; offsets2.setIJ(0,0,0)\r
- offsets2[1:]=DataArrayIdType([int(o) for o in offsets])\r
- offsets3=offsets2.deltaShiftIndex() ; offsets2=offsets3.deepCopy() ; offsets3+=1 ; offsets3.computeOffsetsFull()\r
- offsets=offsets3\r
- tmp1=DataArrayIdType(len(offsets2),2) ; tmp1[:,0]=1 ; tmp1[:,1]=offsets2 ; tmp1.rearrange(1) ; tmp1.computeOffsetsFull()\r
- tmp1=DataArrayIdType.Range(1,2*len(offsets2),2).buildExplicitArrByRanges(tmp1)\r
- conn=np.memmap(fd,dtype=rd._type_conn,mode='r',offset=ref+rd._off_conn,shape=(connLgth,))\r
- conn=self.__swapIfNecessary(rd._bo,conn)\r
- types=np.array(types,dtype='int32') ; types=DataArrayIdType(types) ; \r
- types.transformWithIndArr(self.VTKTypes_2_MC)\r
- conn2=DataArrayIdType(offsets.back())\r
- conn2[offsets[0:-1]]=types\r
- conn2[tmp1]=DataArrayIdType([int(c) for c in conn])\r
- m.setConnectivity(conn2,offsets,True)\r
- m.checkConsistencyLight() ; mm=MEDFileUMesh() ; mm.setMeshAtLevel(0,m) ; ms.pushMesh(mm)\r
- # Fields on nodes and on cells\r
- for spatialDisc,nbEnt,fields in [(ON_NODES,rd._nb_nodes,rd._node_fields),(ON_CELLS,rd._nb_cells,rd._cell_fields)]: \r
- for name,typ,nbCompo,off in fields:\r
- ff=MEDFileFieldMultiTS()\r
- f=MEDCouplingFieldDouble(spatialDisc,ONE_TIME)\r
- f.setName(name) ; f.setMesh(m)\r
- vals=np.memmap(fd,dtype=typ,mode='r',offset=ref+off,shape=(nbEnt*nbCompo))\r
- vals=self.__swapIfNecessary(rd._bo,vals)\r
- arr=DataArrayDouble(np.array(vals,dtype='float64')) ; arr.rearrange(nbCompo)\r
- f.setArray(arr) ; f.checkConsistencyLight()\r
- f.setTime(self._time[0],self._time[1],0)\r
- ff.appendFieldNoProfileSBT(f)\r
- fs.pushField(ff)\r
- pass\r
- pass\r
- return ret\r
-\r
- def __parseXML(self,fd):\r
- import xml.sax\r
- class VTU_SAX_Reader(xml.sax.ContentHandler):\r
- def __init__(self):\r
- self._loc=None\r
- self._data_array={0:self.DAPoints,1:self.DACells,2:self.DAPointData,3:self.DACellData}\r
- self._node_fields=[]\r
- self._cell_fields=[]\r
- pass\r
- def setLocator(self,loc):\r
- self._loc=loc\r
- def DAPoints(self,attrs):\r
- self._space_dim=int(attrs["NumberOfComponents"])\r
- self._type_coords=str(attrs["type"]).lower()\r
- self._off_coords=int(attrs["offset"])\r
- pass\r
- def DACells(self,attrs):\r
- if attrs["Name"]=="connectivity":\r
- self._type_conn=str(attrs["type"]).lower()\r
- self._off_conn=int(attrs["offset"])\r
- pass\r
- if attrs["Name"]=="offsets":\r
- self._type_off=str(attrs["type"]).lower()\r
- #self._off_off=int(attrs.get("offset", 0))\r
- self._off_off=int(attrs["offset"])\r
- pass\r
- if attrs["Name"]=="types":\r
- self._type_types=str(attrs["type"]).lower()\r
- #self._off_types=int(attrs.get("offset", 0))\r
- self._off_types=int(attrs["offset"])\r
- pass\r
- pass\r
- def DAPointData(self,attrs):\r
- numCompo = int(attrs.get("NumberOfComponents", 1))\r
- offset = int(attrs.get("offset", 0))\r
- self._node_fields.append((str(attrs["Name"]),str(attrs["type"]).lower(),numCompo,offset))\r
- pass\r
- def DACellData(self,attrs):\r
- self._cell_fields.append((str(attrs["Name"]),str(attrs["type"]).lower(),int(attrs["NumberOfComponents"]),int(attrs["offset"])))\r
- pass\r
- def startElement(self,name,attrs):\r
- if name=="VTKFile":\r
- if attrs["type"]!="UnstructuredGrid":\r
- raise Exception("Mismatch between reader (VTU) type and file content !")\r
- self._bo=bool(["LittleEndian","BigEndian"].index(attrs["byte_order"]))\r
- pass\r
- if name=="Piece":\r
- self._nb_cells=int(attrs["NumberOfCells"])\r
- self._nb_nodes=int(attrs["NumberOfPoints"])\r
- return\r
- if name=="Points":\r
- self._tmp=0\r
- return\r
- if name=="Cells":\r
- self._tmp=1\r
- return\r
- if name=="PointData":\r
- self._tmp=2\r
- return\r
- if name=="CellData":\r
- self._tmp=3\r
- return\r
- if name=="DataArray":\r
- self._data_array[self._tmp](attrs)\r
- return\r
- if name=="AppendedData":\r
- if str(attrs["encoding"])=="raw":\r
- raise VTURawReader.NormalException(self._loc.getLineNumber())\r
- else:\r
- print(attrs["encoding"])\r
- raise VTURawReader.NotRawVTUException("The file is not a raw VTU ! Change reader !")\r
- pass\r
- pass\r
- rd=VTU_SAX_Reader()\r
- parser=xml.sax.make_parser()\r
- parser.setContentHandler(rd)\r
- locator=xml.sax.expatreader.ExpatLocator(parser)\r
- rd.setLocator(locator)\r
- isOK=False\r
- try:\r
- parser.parse(fd)\r
- except self.NormalException as e:\r
- isOK=True\r
- fd.seek(0)\r
- for i in range(e.getLineNb()): \r
- l = fd.readline()\r
- ref=fd.tell()+12\r
- pass\r
- if not isOK:\r
- raise Exception("Error in VTURawReader : not a raw format ?")\r
- return ref,rd\r
-\r
- @classmethod\r
- def New(cls,fileName,tim=(0.,0)):\r
- """ Static constructor. """\r
- return VTURawReader(fileName,tim)\r
- pass\r
-\r
- def __init__(self,fileName,tim=(0.,0)):\r
- msg="The time specified in constructor as 2nd arg should be a tuple containing 2 values 1 float and 1 int !"\r
- if not isinstance(tim, tuple):\r
- raise Exception(msg)\r
- if len(tim)!=2:\r
- raise Exception(msg)\r
- if not isinstance(tim[0], float) or not isinstance(tim[1], int):\r
- raise Exception(msg)\r
- self._fileName=fileName\r
- self._time=tim\r
- pass\r
-\r
- def __swapIfNecessary(self,b,arr):\r
- if b:\r
- ret=arr.copy()\r
- ret.byteswap(True)\r
- return ret\r
- else:\r
- return arr\r
- pass\r
- pass\r
