Back to index

salome-paravis  6.5.0
VTKMEDCouplingFieldClient.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2010-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 #include "VTKMEDCouplingFieldClient.hxx"
00021 #include "VTKMEDCouplingMeshClient.hxx"
00022 
00023 #include "vtkErrorCode.h"
00024 #include "vtkCellData.h"
00025 #include "vtkPointData.h"
00026 #include "vtkDoubleArray.h"
00027 #include "vtkUnstructuredGrid.h"
00028 
00029 #include <string>
00030 #include <sstream>
00031 
00032 std::vector<double> ParaMEDMEM2VTK::FillMEDCouplingFieldDoubleInstanceFrom(SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr fieldPtr, vtkDataSet *ret)
00033 {
00034   fieldPtr->Register();
00035   //
00036   SALOME_MED::MEDCouplingMeshCorbaInterface_var meshPtr=fieldPtr->getMesh();
00037   ParaMEDMEM2VTK::FillMEDCouplingMeshInstanceFrom(meshPtr,ret);
00038   meshPtr->UnRegister();
00039   //
00040   std::vector<double> ret2=FillMEDCouplingFieldDoublePartOnly(fieldPtr,ret);
00041   fieldPtr->UnRegister();
00042   //
00043   return ret2;
00044 }
00045 
00046 std::vector<double> ParaMEDMEM2VTK::FillMEDCouplingFieldDoublePartOnly(SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr fieldPtr, vtkDataSet *ret)
00047 {
00048   std::vector<double> ret2;
00049   //
00050   SALOME_TYPES::ListOfLong *tinyL;
00051   SALOME_TYPES::ListOfDouble *tinyD;
00052   SALOME_TYPES::ListOfString *tinyS;
00053   fieldPtr->getTinyInfo(tinyL,tinyD,tinyS);
00054   //
00055   int typeOfField=(*tinyL)[0];// 0:ON_CELLS, 1:ON_NODES, 2:ON_GAUSS_PT, 3:ON_GAUSS_NE
00056   int timeDiscr=(*tinyL)[1];
00057   int nbOfTuples=(*tinyL)[3];
00058   int nbOfComponents=(*tinyL)[4];
00059   std::vector<std::string> comps(nbOfComponents);
00060   for(int i=0;i<nbOfComponents;i++)
00061     {
00062       std::string comp((*tinyS)[i]);
00063       if(!comp.empty())
00064         comps[i]=comp;
00065       else
00066         {
00067           std::ostringstream oss; oss << "Comp" << i;
00068           comps[i]=oss.str();
00069         }
00070     }
00071   std::string name;
00072   if(timeDiscr!=7)
00073     name=((*tinyS)[nbOfComponents]);
00074   else
00075     name=((*tinyS)[2*nbOfComponents]);
00076   //
00077   switch(timeDiscr)
00078     {
00079     case 4://NO_TIME
00080       ret2.resize(1);
00081       ret2[0]=-1;
00082       break;
00083     case 5://ONE_TIME
00084       ret2.resize(1);
00085       ret2[0]=(*tinyD)[1];
00086       break;
00087     case 6://LINEAR_TIME
00088     case 7://CONST_ON_TIME_INTERVAL
00089       ret2.resize(1);
00090       ret2[0]=((*tinyD)[1]+(*tinyD)[2])/2.;
00091       break;
00092     default:
00093       vtkErrorWithObjectMacro(ret,"Unrecognized time discretization of Field ! Not implemented yet !");
00094     }
00095   //
00096   delete tinyL;
00097   delete tinyD;
00098   delete tinyS;
00099   //
00100   vtkDoubleArray *var0=vtkDoubleArray::New();
00101   vtkDoubleArray *var1=0;
00102   double *var0Ptr=0;
00103   double *var1Ptr=0;
00104   var0->SetName(name.c_str());
00105   var0->SetNumberOfComponents(nbOfComponents);
00106   var0->SetNumberOfTuples(nbOfTuples);
00107   for(int i=0;i<nbOfComponents;i++)
00108     var0->SetComponentName(i,comps[i].c_str());
00109   var0Ptr=var0->GetPointer(0);
00110   if(timeDiscr==7)
00111     {
00112       var1->SetNumberOfComponents(nbOfComponents);
00113       var1->SetNumberOfTuples(nbOfTuples);
00114       var1Ptr=var1->GetPointer(0);
00115       std::ostringstream oss; oss << name << "_end_array";
00116       var1->SetName(oss.str().c_str());
00117     }
00118   //
00119   SALOME_TYPES::ListOfLong *bigDataI;
00120   SALOME_TYPES::ListOfDouble2 *bigArr;
00121   fieldPtr->getSerialisationData(bigDataI,bigArr);
00122   delete bigDataI;//for the moment gauss points are not dealt
00123   SALOME_TYPES::ListOfDouble& oneArr=(*bigArr)[0];
00124   int nbVals=oneArr.length();
00125   for(int i=0;i<nbVals;i++)
00126     var0Ptr[i]=oneArr[i];
00127   if(timeDiscr==7)
00128     {
00129       SALOME_TYPES::ListOfDouble& scndArr=(*bigArr)[1];
00130       for(int i=0;i<nbVals;i++)
00131         var1Ptr[i]=scndArr[i];
00132     }
00133   delete bigArr;
00134   //
00135   switch(typeOfField)
00136     {
00137     case 0://ON_CELLS
00138       ret->GetCellData()->AddArray(var0);
00139       if(timeDiscr==7)
00140         ret->GetCellData()->AddArray(var1);
00141       break;
00142     case 1://ON_NODES
00143       ret->GetPointData()->AddArray(var0);
00144       if(timeDiscr==7)
00145         ret->GetCellData()->AddArray(var1);
00146       break;
00147     default:
00148       vtkErrorWithObjectMacro(ret,"Not implemented yet for gauss fields and gauss on elements fields !");
00149     }
00150   var0->Delete();
00151   if(timeDiscr==7)
00152     var1->Delete();
00153   //
00154   return ret2;
00155 }
00156 
00157 vtkDataSet *ParaMEDMEM2VTK::BuildFullyFilledFromMEDCouplingFieldDoubleInstance(SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr fieldPtr,
00158                                                                                std::vector<double>& times)
00159 {
00160   fieldPtr->Register();
00161   SALOME_MED::MEDCouplingMeshCorbaInterface_var meshPtr=fieldPtr->getMesh();
00162   bool dummy;//VTK bug
00163   vtkDataSet *ret=ParaMEDMEM2VTK::BuildFromMEDCouplingMeshInstance(meshPtr,dummy);//VTK bug
00164   meshPtr->UnRegister();
00165   //
00166   std::vector<double> ret2=FillMEDCouplingFieldDoublePartOnly(fieldPtr,ret);
00167   times=ret2;
00168   //
00169   fieldPtr->UnRegister();
00170   return ret;
00171 }
00172 
00173 vtkDoubleArray *ParaMEDMEM2VTK::BuildFromMEDCouplingFieldDoubleArr(SALOME_MED::DataArrayDoubleCorbaInterface_ptr dadPtr)
00174 {
00175   vtkDoubleArray *ret=vtkDoubleArray::New();
00176   //
00177   SALOME_TYPES::ListOfLong *tinyL=0;
00178   SALOME_TYPES::ListOfDouble *bigD=0;
00179   SALOME_TYPES::ListOfString *tinyS=0;
00180   //
00181   dadPtr->getTinyInfo(tinyL,tinyS);
00182   int nbOfTuples=(*tinyL)[0];
00183   int nbOfCompo=(*tinyL)[1];
00184   std::string name((*tinyS)[0]);
00185   std::vector<std::string> comps(nbOfCompo);
00186   for(int i=0;i<nbOfCompo;i++)
00187     comps[i]=(*tinyS)[i+1];
00188   delete tinyL;
00189   delete tinyS;
00190   //
00191   ret->SetName(name.c_str());
00192   ret->SetNumberOfComponents(nbOfCompo);
00193   ret->SetNumberOfTuples(nbOfTuples);
00194   for(int i=0;i<nbOfCompo;i++)
00195     {
00196       if(!comps[i].empty())
00197         ret->SetComponentName(i,comps[i].c_str());
00198       else
00199         {
00200           std::ostringstream oss; oss << "Comp" << i;
00201           ret->SetComponentName(i,oss.str().c_str());
00202         }
00203     }
00204   int nbElems=nbOfCompo*nbOfTuples;
00205   double *pt=ret->GetPointer(0);
00206   dadPtr->getSerialisationData(bigD);
00207   for(int i=0;i<nbElems;i++)
00208     pt[i]=(*bigD)[i];
00209   delete bigD;
00210   //
00211   return ret;
00212 }