Back to index

salome-paravis  6.5.0
Classes | Functions
ParaMEDMEM2VTK Namespace Reference

Classes

class  TinyInfoOnField
 Stores all info on field without consideration of time. More...
protocol  ParaMEDMEM2VTK_EXPORT

Functions

void FillMEDCouplingCMeshInstanceFrom (SALOME_MED::MEDCouplingCMeshCorbaInterface_ptr meshPtr, vtkRectilinearGrid *ret)
std::vector< double > FillMEDCouplingFieldDoubleInstanceFrom (SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr fieldPtr, vtkDataSet *ret)
std::vector< double > FillMEDCouplingFieldDoublePartOnly (SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr fieldPtr, vtkDataSet *ret)
ParaMEDMEM2VTK_EXPORT vtkDataSet * BuildFullyFilledFromMEDCouplingFieldDoubleInstance (SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr fieldPtr, std::vector< double > &times)
ParaMEDMEM2VTK_EXPORT
vtkDoubleArray * 
BuildFromMEDCouplingFieldDoubleArr (SALOME_MED::DataArrayDoubleCorbaInterface_ptr dadPtr)
void FillMEDCouplingMeshInstanceFrom (SALOME_MED::MEDCouplingMeshCorbaInterface_ptr meshPtr, vtkDataSet *ret)
ParaMEDMEM2VTK_EXPORT vtkDataSet * BuildFromMEDCouplingMeshInstance (SALOME_MED::MEDCouplingMeshCorbaInterface_ptr meshPtr, bool &isPolyh)
void FillMEDCouplingUMeshInstanceFrom (SALOME_MED::MEDCouplingUMeshCorbaInterface_ptr meshPtr, vtkUnstructuredGrid *ret, bool &isPolyh)
ParaMEDMEM2VTK_EXPORT
std::vector< double > 
FillMEDCouplingParaFieldDoubleInstanceFrom (SALOME_MED::ParaMEDCouplingFieldDoubleCorbaInterface_ptr fieldPtr, int begin, int end, vtkMultiBlockDataSet *ret)

Class Documentation

class ParaMEDMEM2VTK::TinyInfoOnField

Stores all info on field without consideration of time.

Definition at line 37 of file VTKMEDCouplingMultiFieldsClient.hxx.

Collaboration diagram for ParaMEDMEM2VTK::TinyInfoOnField:
Class Members
string _name
int _type

Function Documentation

vtkDoubleArray * ParaMEDMEM2VTK::BuildFromMEDCouplingFieldDoubleArr ( SALOME_MED::DataArrayDoubleCorbaInterface_ptr  dadPtr)

Definition at line 173 of file VTKMEDCouplingFieldClient.cxx.

{
  vtkDoubleArray *ret=vtkDoubleArray::New();
  //
  SALOME_TYPES::ListOfLong *tinyL=0;
  SALOME_TYPES::ListOfDouble *bigD=0;
  SALOME_TYPES::ListOfString *tinyS=0;
  //
  dadPtr->getTinyInfo(tinyL,tinyS);
  int nbOfTuples=(*tinyL)[0];
  int nbOfCompo=(*tinyL)[1];
  std::string name((*tinyS)[0]);
  std::vector<std::string> comps(nbOfCompo);
  for(int i=0;i<nbOfCompo;i++)
    comps[i]=(*tinyS)[i+1];
  delete tinyL;
  delete tinyS;
  //
  ret->SetName(name.c_str());
  ret->SetNumberOfComponents(nbOfCompo);
  ret->SetNumberOfTuples(nbOfTuples);
  for(int i=0;i<nbOfCompo;i++)
    {
      if(!comps[i].empty())
        ret->SetComponentName(i,comps[i].c_str());
      else
        {
          std::ostringstream oss; oss << "Comp" << i;
          ret->SetComponentName(i,oss.str().c_str());
        }
    }
  int nbElems=nbOfCompo*nbOfTuples;
  double *pt=ret->GetPointer(0);
  dadPtr->getSerialisationData(bigD);
  for(int i=0;i<nbElems;i++)
    pt[i]=(*bigD)[i];
  delete bigD;
  //
  return ret;
}
vtkDataSet * ParaMEDMEM2VTK::BuildFromMEDCouplingMeshInstance ( SALOME_MED::MEDCouplingMeshCorbaInterface_ptr  meshPtr,
bool &  isPolyh 
)

Definition at line 62 of file VTKMEDCouplingMeshClient.cxx.

{
  SALOME_MED::MEDCouplingUMeshCorbaInterface_var umeshPtr=SALOME_MED::MEDCouplingUMeshCorbaInterface::_narrow(meshPtr);
  if(!CORBA::is_nil(umeshPtr))
    {
      vtkUnstructuredGrid *ret1=vtkUnstructuredGrid::New();
      ParaMEDMEM2VTK::FillMEDCouplingUMeshInstanceFrom(umeshPtr,ret1,isPolyh);
      return ret1;
    }
  SALOME_MED::MEDCouplingCMeshCorbaInterface_var cmeshPtr=SALOME_MED::MEDCouplingCMeshCorbaInterface::_narrow(meshPtr);
  if(!CORBA::is_nil(cmeshPtr))
    {
      vtkRectilinearGrid *ret1=vtkRectilinearGrid::New();
      ParaMEDMEM2VTK::FillMEDCouplingCMeshInstanceFrom(cmeshPtr,ret1);
      return ret1;
    }
  vtkOutputWindowDisplayErrorText("Error : CORBA mesh type ! Mesh type not managed #2 !");
  return 0;
}

Here is the call graph for this function:

Here is the caller graph for this function:

vtkDataSet * ParaMEDMEM2VTK::BuildFullyFilledFromMEDCouplingFieldDoubleInstance ( SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr  fieldPtr,
std::vector< double > &  times 
)

Definition at line 157 of file VTKMEDCouplingFieldClient.cxx.

{
  fieldPtr->Register();
  SALOME_MED::MEDCouplingMeshCorbaInterface_var meshPtr=fieldPtr->getMesh();
  bool dummy;//VTK bug
  vtkDataSet *ret=ParaMEDMEM2VTK::BuildFromMEDCouplingMeshInstance(meshPtr,dummy);//VTK bug
  meshPtr->UnRegister();
  //
  std::vector<double> ret2=FillMEDCouplingFieldDoublePartOnly(fieldPtr,ret);
  times=ret2;
  //
  fieldPtr->UnRegister();
  return ret;
}

Here is the call graph for this function:

Here is the caller graph for this function:

void ParaMEDMEM2VTK::FillMEDCouplingCMeshInstanceFrom ( SALOME_MED::MEDCouplingCMeshCorbaInterface_ptr  meshPtr,
vtkRectilinearGrid *  ret 
)

Definition at line 29 of file VTKMEDCouplingCMeshClient.cxx.

{
  meshPtr->Register();
  SALOME_TYPES::ListOfDouble *tinyD;
  SALOME_TYPES::ListOfLong *tinyI;
  SALOME_TYPES::ListOfString *tinyS;
  meshPtr->getTinyInfo(tinyD,tinyI,tinyS);
  int sizePerAxe[3];
  sizePerAxe[0]=(*tinyI)[0];
  sizePerAxe[1]=(*tinyI)[1];
  sizePerAxe[2]=(*tinyI)[2];
  ret->SetDimensions(sizePerAxe[0],sizePerAxe[1],sizePerAxe[2]);
  delete tinyI;
  delete tinyS;
  SALOME_TYPES::ListOfDouble *bigD;
  meshPtr->getSerialisationData(tinyI,bigD);
  delete tinyI;
  int offset=0;
  vtkDoubleArray *da=0;
  if(sizePerAxe[0]>0)
    {
      da=vtkDoubleArray::New();
      da->SetNumberOfTuples(sizePerAxe[0]);
      da->SetNumberOfComponents(1);
      double *pt=da->GetPointer(0);
      for(int i=0;i<sizePerAxe[0];i++)
        pt[i]=(*bigD)[i];
      ret->SetXCoordinates(da);
      da->Delete();
      offset+=sizePerAxe[0];
    }
  if(sizePerAxe[1]>0)
    {
      da=vtkDoubleArray::New();
      da->SetNumberOfTuples(sizePerAxe[1]);
      da->SetNumberOfComponents(1);
      double *pt=da->GetPointer(0);
      for(int i=0;i<sizePerAxe[1];i++)
        pt[i]=(*bigD)[offset+i];
      ret->SetYCoordinates(da);
      da->Delete();
      offset+=sizePerAxe[1];
    }
  if(sizePerAxe[2]>0)
    {
      da=vtkDoubleArray::New();
      da->SetNumberOfTuples(sizePerAxe[2]);
      da->SetNumberOfComponents(1);
      double *pt=da->GetPointer(0);
      for(int i=0;i<sizePerAxe[2];i++)
        pt[i]=(*bigD)[offset+i];
      ret->SetZCoordinates(da);
      da->Delete();
    }
  delete bigD;
  meshPtr->UnRegister();
}

Here is the caller graph for this function:

std::vector< double > ParaMEDMEM2VTK::FillMEDCouplingFieldDoubleInstanceFrom ( SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr  fieldPtr,
vtkDataSet *  ret 
)

Definition at line 32 of file VTKMEDCouplingFieldClient.cxx.

{
  fieldPtr->Register();
  //
  SALOME_MED::MEDCouplingMeshCorbaInterface_var meshPtr=fieldPtr->getMesh();
  ParaMEDMEM2VTK::FillMEDCouplingMeshInstanceFrom(meshPtr,ret);
  meshPtr->UnRegister();
  //
  std::vector<double> ret2=FillMEDCouplingFieldDoublePartOnly(fieldPtr,ret);
  fieldPtr->UnRegister();
  //
  return ret2;
}

Here is the call graph for this function:

std::vector< double > ParaMEDMEM2VTK::FillMEDCouplingFieldDoublePartOnly ( SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_ptr  fieldPtr,
vtkDataSet *  ret 
)

Definition at line 46 of file VTKMEDCouplingFieldClient.cxx.

{
  std::vector<double> ret2;
  //
  SALOME_TYPES::ListOfLong *tinyL;
  SALOME_TYPES::ListOfDouble *tinyD;
  SALOME_TYPES::ListOfString *tinyS;
  fieldPtr->getTinyInfo(tinyL,tinyD,tinyS);
  //
  int typeOfField=(*tinyL)[0];// 0:ON_CELLS, 1:ON_NODES, 2:ON_GAUSS_PT, 3:ON_GAUSS_NE
  int timeDiscr=(*tinyL)[1];
  int nbOfTuples=(*tinyL)[3];
  int nbOfComponents=(*tinyL)[4];
  std::vector<std::string> comps(nbOfComponents);
  for(int i=0;i<nbOfComponents;i++)
    {
      std::string comp((*tinyS)[i]);
      if(!comp.empty())
        comps[i]=comp;
      else
        {
          std::ostringstream oss; oss << "Comp" << i;
          comps[i]=oss.str();
        }
    }
  std::string name;
  if(timeDiscr!=7)
    name=((*tinyS)[nbOfComponents]);
  else
    name=((*tinyS)[2*nbOfComponents]);
  //
  switch(timeDiscr)
    {
    case 4://NO_TIME
      ret2.resize(1);
      ret2[0]=-1;
      break;
    case 5://ONE_TIME
      ret2.resize(1);
      ret2[0]=(*tinyD)[1];
      break;
    case 6://LINEAR_TIME
    case 7://CONST_ON_TIME_INTERVAL
      ret2.resize(1);
      ret2[0]=((*tinyD)[1]+(*tinyD)[2])/2.;
      break;
    default:
      vtkErrorWithObjectMacro(ret,"Unrecognized time discretization of Field ! Not implemented yet !");
    }
  //
  delete tinyL;
  delete tinyD;
  delete tinyS;
  //
  vtkDoubleArray *var0=vtkDoubleArray::New();
  vtkDoubleArray *var1=0;
  double *var0Ptr=0;
  double *var1Ptr=0;
  var0->SetName(name.c_str());
  var0->SetNumberOfComponents(nbOfComponents);
  var0->SetNumberOfTuples(nbOfTuples);
  for(int i=0;i<nbOfComponents;i++)
    var0->SetComponentName(i,comps[i].c_str());
  var0Ptr=var0->GetPointer(0);
  if(timeDiscr==7)
    {
      var1->SetNumberOfComponents(nbOfComponents);
      var1->SetNumberOfTuples(nbOfTuples);
      var1Ptr=var1->GetPointer(0);
      std::ostringstream oss; oss << name << "_end_array";
      var1->SetName(oss.str().c_str());
    }
  //
  SALOME_TYPES::ListOfLong *bigDataI;
  SALOME_TYPES::ListOfDouble2 *bigArr;
  fieldPtr->getSerialisationData(bigDataI,bigArr);
  delete bigDataI;//for the moment gauss points are not dealt
  SALOME_TYPES::ListOfDouble& oneArr=(*bigArr)[0];
  int nbVals=oneArr.length();
  for(int i=0;i<nbVals;i++)
    var0Ptr[i]=oneArr[i];
  if(timeDiscr==7)
    {
      SALOME_TYPES::ListOfDouble& scndArr=(*bigArr)[1];
      for(int i=0;i<nbVals;i++)
        var1Ptr[i]=scndArr[i];
    }
  delete bigArr;
  //
  switch(typeOfField)
    {
    case 0://ON_CELLS
      ret->GetCellData()->AddArray(var0);
      if(timeDiscr==7)
        ret->GetCellData()->AddArray(var1);
      break;
    case 1://ON_NODES
      ret->GetPointData()->AddArray(var0);
      if(timeDiscr==7)
        ret->GetCellData()->AddArray(var1);
      break;
    default:
      vtkErrorWithObjectMacro(ret,"Not implemented yet for gauss fields and gauss on elements fields !");
    }
  var0->Delete();
  if(timeDiscr==7)
    var1->Delete();
  //
  return ret2;
}

Here is the caller graph for this function:

void ParaMEDMEM2VTK::FillMEDCouplingMeshInstanceFrom ( SALOME_MED::MEDCouplingMeshCorbaInterface_ptr  meshPtr,
vtkDataSet *  ret 
)

Definition at line 32 of file VTKMEDCouplingMeshClient.cxx.

{
  SALOME_MED::MEDCouplingUMeshCorbaInterface_var umeshPtr=SALOME_MED::MEDCouplingUMeshCorbaInterface::_narrow(meshPtr);
  if(!CORBA::is_nil(umeshPtr))
    {
      vtkUnstructuredGrid *ret1=vtkUnstructuredGrid::SafeDownCast(ret);
      if(!ret1)
        {
          vtkErrorWithObjectMacro(ret,"Internal error in ParaMEDCorba plugin : mismatch between VTK type and CORBA type UMesh !");
          return ;
        }
      bool dummy;//VTK bug
      ParaMEDMEM2VTK::FillMEDCouplingUMeshInstanceFrom(umeshPtr,ret1,dummy);//VTK bug
      return ;
    }
  SALOME_MED::MEDCouplingCMeshCorbaInterface_var cmeshPtr=SALOME_MED::MEDCouplingCMeshCorbaInterface::_narrow(meshPtr);
  if(!CORBA::is_nil(cmeshPtr))
    {
      vtkRectilinearGrid *ret1=vtkRectilinearGrid::SafeDownCast(ret);
      if(!ret1)
        {
          vtkErrorWithObjectMacro(ret,"Internal error in ParaMEDCorba plugin : mismatch between VTK type and CORBA type CMesh !");
          return ;
        }
      ParaMEDMEM2VTK::FillMEDCouplingCMeshInstanceFrom(cmeshPtr,ret1);
      return ;
    }
  vtkErrorWithObjectMacro(ret,"Error : CORBA mesh type ! Mesh type not managed !");
}

Here is the call graph for this function:

Here is the caller graph for this function:

std::vector< double > ParaMEDMEM2VTK::FillMEDCouplingParaFieldDoubleInstanceFrom ( SALOME_MED::ParaMEDCouplingFieldDoubleCorbaInterface_ptr  fieldPtr,
int  begin,
int  end,
vtkMultiBlockDataSet *  ret 
)

Definition at line 29 of file VTKParaMEDFieldClient.cxx.

{
  std::vector<double> ret2;
  int nbOfParts=end-begin;
  Engines::IORTab *allSlices=fieldPtr->tior();
  vtkMultiBlockDataGroupFilter *tmp=vtkMultiBlockDataGroupFilter::New();
  for(int i=0;i<nbOfParts;i++)
    {
      CORBA::Object_ptr obj=(*allSlices)[i+begin];
      SALOME_MED::MEDCouplingFieldDoubleCorbaInterface_var fieldCorba=SALOME_MED::MEDCouplingFieldDoubleCorbaInterface::_narrow(obj);
      std::vector<double> times;
      vtkDataSet *part=ParaMEDMEM2VTK::BuildFullyFilledFromMEDCouplingFieldDoubleInstance(fieldCorba,times);
      tmp->AddInputConnection(part->GetProducerPort());
      part->Delete();
    }
  vtkCompositeDataToUnstructuredGridFilter *tmp2=vtkCompositeDataToUnstructuredGridFilter::New();
  tmp2->AddInputConnection(tmp->GetOutput()->GetProducerPort());
  tmp2->Update();
  //
  vtkUnstructuredGrid *ret3=tmp2->GetOutput();
  ret->SetBlock(0,ret3);
  //
  tmp->Delete();
  tmp2->Delete();
  delete allSlices;
  return ret2;
}

Here is the call graph for this function:

Here is the caller graph for this function:

void ParaMEDMEM2VTK::FillMEDCouplingUMeshInstanceFrom ( SALOME_MED::MEDCouplingUMeshCorbaInterface_ptr  meshPtr,
vtkUnstructuredGrid *  ret,
bool &  isPolyh 
)

Definition at line 35 of file VTKMEDCouplingUMeshClient.cxx.

{
  meshPtr->Register();
  //
  SALOME_TYPES::ListOfDouble *tinyD;
  SALOME_TYPES::ListOfLong *tinyI;
  SALOME_TYPES::ListOfString *tinyS;
  meshPtr->getTinyInfo(tinyD,tinyI,tinyS);
  //
  int spaceDim=(*tinyI)[1];
  int nbOfNodes=(*tinyI)[2];
  int meshDim=(*tinyI)[5];
  int nbOfCells=(*tinyI)[6];
  int meshLength=(*tinyI)[7];
  std::string name((*tinyS)[0]);
  //std::vector<std::string> compoNames(spaceDim);
  //for(int i=0;i<spaceDim;i++)
  //  compoNames[i]=(*tinyS)[i+4];
  delete tinyD;
  delete tinyI;
  delete tinyS;
  //
  ret->Initialize();
  ret->Allocate(nbOfCells);
  vtkPoints *points=vtkPoints::New();
  vtkDoubleArray *da=vtkDoubleArray::New();
  da->SetNumberOfComponents(3);
  da->SetNumberOfTuples(nbOfNodes);
  double *pts=da->GetPointer(0);
  //
  SALOME_TYPES::ListOfLong *a1Corba;
  SALOME_TYPES::ListOfDouble *a2Corba;
  meshPtr->getSerialisationData(a1Corba,a2Corba);
  if(spaceDim==3)
    {
      int myLgth=a2Corba->length();
      for(int i=0;i<myLgth;i++)
        *pts++=(*a2Corba)[i];
    }
  else
    {
      int offset=3-spaceDim;
      for(int i=0;i<nbOfNodes;i++)
        {
          for(int j=0;j<spaceDim;j++)
            *pts++=(*a2Corba)[spaceDim*i+j];
          std::fill(pts,pts+offset,0.);
          pts+=offset;
        }
    }
  //
  vtkIdType *tmp=new vtkIdType[1000];
  isPolyh=false;
  for(int i=0;i<nbOfCells;i++)
    {
      int pos=(*a1Corba)[i];
      int pos1=(*a1Corba)[i+1];
      int nbOfNodeInCurCell=pos1-pos-1;
      int typeOfCell=(*a1Corba)[pos+nbOfCells+1];
      for(int j=0;j<nbOfNodeInCurCell;j++)
        tmp[j]=(*a1Corba)[pos+1+j+nbOfCells+1];
      int vtkType=ParaMEDMEM2VTKTypeTraducer[typeOfCell];
      if(vtkType!=42)
        ret->InsertNextCell(vtkType,nbOfNodeInCurCell,tmp);
      else
        {//polyhedron
          isPolyh=true;
          std::set<vtkIdType> s(tmp,tmp+nbOfNodeInCurCell);
          vtkSmartPointer<vtkCellArray> faces=vtkSmartPointer<vtkCellArray>::New();
          int nbOfFaces=std::count(tmp,tmp+nbOfNodeInCurCell,-1)+1;
          vtkIdType *work=tmp;
          for(int i=0;i<nbOfFaces;i++)
            {
              vtkIdType *work2=std::find(work,tmp+nbOfNodeInCurCell,-1);
              int nbOfNodesInFace=std::distance(work,work2);
              faces->InsertNextCell(nbOfNodesInFace,work);
              work=work2+1;
            }
          s.erase(-1);
          std::vector<vtkIdType> v(s.begin(),s.end());
          ret->InsertNextCell(VTK_POLYHEDRON,v.size(),&v[0],nbOfFaces,faces->GetPointer());
        }
    }
  delete [] tmp;
  //
  delete a1Corba;
  delete a2Corba;
  //
  ret->SetPoints(points);
  points->Delete();
  points->SetData(da);
  da->Delete();
  //
  meshPtr->UnRegister();
}

Here is the caller graph for this function: