Back to index

salome-paravis  6.5.0
VTKMEDCouplingMultiFieldsClient.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 // To access to vtkUnstructuredGrid::Faces and FaceLocations
00021 #define protected public
00022 
00023 #include "VTKMEDCouplingMultiFieldsClient.hxx"
00024 #include "VTKMEDCouplingMeshClient.hxx"
00025 #include "VTKMEDCouplingFieldClient.hxx"
00026 
00027 #include "vtkUnstructuredGrid.h"
00028 #include "vtkRectilinearGrid.h"
00029 #include "vtkDoubleArray.h"
00030 #include "vtkErrorCode.h"
00031 #include "vtkCellData.h"
00032 #include "vtkIdTypeArray.h"
00033 #include "vtkPointData.h"
00034 
00035 #include <sstream>
00036 #include <iterator>
00037 #include <algorithm>
00038 #include <functional>
00039 
00040 const double ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::EPS_TIME=1e-7;
00041 
00042 ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::MEDCouplingMultiFieldsFetcher(int bufferingPolicy,
00043                                                                     SALOME_MED::MEDCouplingMultiFieldsCorbaInterface_ptr mfieldsPtr):_effective_pol(bufferingPolicy),_mfields_ptr_released(false)
00044 {
00045   _mfields_ptr=SALOME_MED::MEDCouplingMultiFieldsCorbaInterface::_duplicate(mfieldsPtr);
00046   _mfields_ptr->Register();
00047 }
00048 
00049 ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::~MEDCouplingMultiFieldsFetcher()
00050 {
00051   for(std::vector<vtkDataSet *>::iterator it=_meshes.begin();it!=_meshes.end();it++)
00052     {
00053       if(*it)
00054        (*it)->Delete();
00055     }
00056   for(std::vector<vtkDoubleArray *>::iterator it2=_arrays.begin();it2!=_arrays.end();it2++)
00057     {
00058       if(*it2)
00059        (*it2)->Delete();
00060     }
00061   if(!_mfields_ptr_released)
00062     _mfields_ptr->UnRegister();
00063 }
00064 
00065 std::vector<double> ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::getTimeStepsForPV()
00066 {
00067   retrievesMainTinyInfo();
00068   int nbOfFields=_mesh_id_per_field.size();
00069   //
00070   _time_label_per_field.resize(nbOfFields);
00071   SALOME_MED::MEDCouplingFieldOverTimeCorbaInterface_var fotPtr=SALOME_MED::MEDCouplingFieldOverTimeCorbaInterface::_narrow(_mfields_ptr);
00072   if(CORBA::is_nil(fotPtr))
00073     {
00074       for(CORBA::Long i=0;i<nbOfFields;i++)
00075        _time_label_per_field[i]=(double)i;
00076     }
00077   else
00078     {
00079       double tmp=0.;
00080       for(CORBA::Long i=0;i<nbOfFields;i++)
00081        {
00082          if(!_time_def_per_field[i].empty())
00083            _time_label_per_field[i]=_time_def_per_field[i].front();
00084          else
00085            _time_label_per_field[i]=tmp++;
00086        }
00087     }
00088   return _time_label_per_field;
00089 }
00090 
00091 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::fetchRegardingPolicy()
00092 {
00093   if(_effective_pol>=10)
00094     {
00095       fetchAll();
00096       return ;
00097     }
00098   if(_effective_pol>=1 && _effective_pol<=9)
00099     {
00100       fetchMeshes();
00101       return ;
00102     }
00103 }
00104 
00105 vtkDataSet *ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::buildDataSetOnTime(double time)
00106 {
00107   int fieldId=getPosGivenTimeLabel(time);
00108   if(fieldId<0)
00109     return 0;
00110   fetchDataIfNeeded(fieldId);
00111   int meshId=_mesh_id_per_field[fieldId];
00112   vtkDataSet *ret0=_meshes[meshId];
00113   std::string clsName=ret0->GetClassName();
00114   if(clsName=="vtkUnstructuredGrid")
00115     {
00116       vtkUnstructuredGrid *ret1=vtkUnstructuredGrid::New();
00117       ret1->DeepCopy(ret0);
00118       if(_is_meshes_polyhedron[meshId])//bug VTK polyhedron
00119         {//bug VTK polyhedron part
00120           ret1->Faces->UnRegister(ret1);
00121           ret1->Faces=vtkIdTypeArray::New();
00122           ret1->Faces->DeepCopy(((vtkUnstructuredGrid *)ret0)->GetFaces());
00123           ret1->Faces->Register(ret1);
00124           ret1->Faces->Delete();
00125           ret1->FaceLocations->UnRegister(ret1);
00126           ret1->FaceLocations=vtkIdTypeArray::New();
00127           ret1->FaceLocations->DeepCopy(((vtkUnstructuredGrid *)ret0)->GetFaceLocations());
00128           ret1->FaceLocations->Register(ret1);
00129           ret1->FaceLocations->Delete();
00130         }//end bug VTK polyhedron part
00131       appendFieldValueOnAlreadyFetchedData(ret1,fieldId);
00132       applyBufferingPolicy();
00133       return ret1;
00134     }
00135   if(clsName=="vtkRectilinearGrid")
00136     {
00137       vtkRectilinearGrid *ret1=vtkRectilinearGrid::New();
00138       ret1->DeepCopy(ret0);
00139       appendFieldValueOnAlreadyFetchedData(ret1,fieldId);
00140       applyBufferingPolicy();
00141       return ret1;
00142     }
00143   return 0;
00144 }
00145 
00146 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::retrievesMainTinyInfo()
00147 {
00148   SALOME_TYPES::ListOfLong *tinyL=0;
00149   SALOME_TYPES::ListOfDouble *tinyD=0;
00150   SALOME_TYPES::ListOfString *tinyS=0;
00151   //
00152   CORBA::Long nbOfArrays;
00153   CORBA::Long nbOfFields;
00154   CORBA::Long nbOfMeshes=_mfields_ptr->getMainTinyInfo(tinyL,tinyD,nbOfArrays,nbOfFields);
00155   int sz=(*tinyL)[0];//nbOfFields
00156   //int sz2=(*tinyL)[1];//sigma(nbOfArraysPerField)
00157   _time_discr_per_field.resize(sz);//4 : NO_TIME  5:ONE_TIME 6:LINEAR_TIME 7:CONST_ON_TIME_INTERVAL
00158   _mesh_id_per_field.resize(sz);
00159   _array_ids_per_field.resize(sz);
00160   _time_def_per_field.resize(sz);
00161   int offsetTime=0;
00162   int offsetArrays=0;
00163   for(int i=0;i<sz;i++)
00164     {
00165       _mesh_id_per_field[i]=(*tinyL)[3+i];
00166       int nbOfArrayForCurField=(*tinyL)[sz+3+i];
00167       _array_ids_per_field[i].resize(nbOfArrayForCurField);
00168       for(int k=0;k<nbOfArrayForCurField;k++)
00169        _array_ids_per_field[i][k]=(*tinyL)[5*sz+3+offsetArrays+k];
00170       _time_discr_per_field[i]=(*tinyL)[2*sz+3+i];
00171       int nbOfTimeSpot=(*tinyL)[3*sz+3+i]-1;//-1 because time precision is not useful here.
00172       _time_def_per_field[i].resize(nbOfTimeSpot);
00173       for(int j=0;j<nbOfTimeSpot;j++)
00174        _time_def_per_field[i][j]=(*tinyD)[offsetTime+1+j];
00175       offsetTime+=nbOfTimeSpot+1;
00176       offsetArrays+=nbOfArrayForCurField;
00177     }
00178   delete tinyL;
00179   delete tinyD;
00180   //
00181   _meshes.resize(nbOfMeshes+1);
00182   _is_meshes_polyhedron.resize(nbOfMeshes+1);
00183   _arrays.resize(nbOfArrays+1);
00184   //
00185   _info_per_field.resize(nbOfFields);
00186   for(int i=0;i<nbOfFields;i++)
00187     {
00188       _mfields_ptr->getTinyInfo(i,tinyL,tinyD,tinyS);
00189       _info_per_field[i]._type=(*tinyL)[0];
00190       _info_per_field[i]._name=(*tinyS)[0];
00191       delete tinyL;
00192       delete tinyD;
00193       delete tinyS;
00194     }
00195 }
00196 
00197 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::fetchAll()
00198 {
00199   fetchMeshes();
00200   int nbOfArrays=_arrays.size();
00201   for(int i=0;i<nbOfArrays;i++)
00202     {
00203       SALOME_MED::DataArrayDoubleCorbaInterface_var daPtr=_mfields_ptr->getArray(i);
00204       if(_arrays[i])
00205        _arrays[i]->Delete();
00206       _arrays[i]=ParaMEDMEM2VTK::BuildFromMEDCouplingFieldDoubleArr(daPtr);
00207       daPtr->UnRegister();
00208     }
00209   unregisterRemoteServantIfAllFetched();
00210 }
00211 
00215 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::fetchMeshes()
00216 {
00217   int nbOfMeshes=_meshes.size();
00218   for(int i=0;i<nbOfMeshes;i++)
00219     {
00220       SALOME_MED::MEDCouplingMeshCorbaInterface_var mPtr=_mfields_ptr->getMeshWithId(i);
00221       if(_meshes[i])
00222        _meshes[i]->Delete();
00223       bool polyh=false;//bug VTK
00224       _meshes[i]=ParaMEDMEM2VTK::BuildFromMEDCouplingMeshInstance(mPtr,polyh);//bug VTK
00225       _is_meshes_polyhedron[i]=polyh;//bug VTK
00226       mPtr->UnRegister();
00227     }
00228   unregisterRemoteServantIfAllFetched();
00229 }
00230 
00235 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::fetchDataIfNeeded(int fieldId)
00236 {
00237   std::vector<int> arrayIds=_array_ids_per_field[fieldId];
00238   int meshId=_mesh_id_per_field[fieldId];
00239   if(!_meshes[meshId])
00240     {
00241       SALOME_MED::MEDCouplingMeshCorbaInterface_var mPtr=_mfields_ptr->getMeshWithId(meshId);
00242       bool polyh=false;//bug VTK
00243       _meshes[meshId]=ParaMEDMEM2VTK::BuildFromMEDCouplingMeshInstance(mPtr,polyh);//bug VTK
00244       _is_meshes_polyhedron[meshId]=polyh;//bug VTK
00245       mPtr->UnRegister();
00246     }
00247   for(std::vector<int>::const_iterator it=arrayIds.begin();it!=arrayIds.end();it++)
00248     {
00249       if(!_arrays[*it])
00250         {
00251           SALOME_MED::DataArrayDoubleCorbaInterface_var daPtr=_mfields_ptr->getArray(*it);
00252           _arrays[*it]=ParaMEDMEM2VTK::BuildFromMEDCouplingFieldDoubleArr(daPtr);
00253           daPtr->UnRegister();
00254         }
00255     }
00256   unregisterRemoteServantIfAllFetched();
00257 }
00258 
00259 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::unregisterRemoteServantIfAllFetched()
00260 {
00261   for(std::vector<vtkDataSet *>::iterator it=_meshes.begin();it!=_meshes.end();it++)
00262     {
00263       if((*it)==0)
00264        return ;
00265     }
00266   for(std::vector<vtkDoubleArray *>::iterator it2=_arrays.begin();it2!=_arrays.end();it2++)
00267     {
00268       if((*it2)==0)
00269        return ;
00270     }
00271   if(!_mfields_ptr_released)
00272     {
00273       _mfields_ptr_released=true;
00274       _mfields_ptr->UnRegister();
00275     }
00276 }
00277 
00278 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::applyBufferingPolicy()
00279 {
00280   if(_effective_pol==0)
00281     {//
00282       for(std::vector<vtkDataSet *>::iterator it=_meshes.begin();it!=_meshes.end();it++)
00283         {
00284           if(*it)
00285             {
00286               (*it)->Delete();
00287               *it=0;
00288             }
00289         }
00290       for(std::vector<vtkDoubleArray *>::iterator it2=_arrays.begin();it2!=_arrays.end();it2++)
00291         {
00292           if(*it2)
00293             {
00294               (*it2)->Delete();
00295               *it2=0;
00296             }
00297         }
00298     }
00299   //else nothing to do let the plugin bufferize
00300 }
00301 
00302 void ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::appendFieldValueOnAlreadyFetchedData(vtkDataSet *ds, int fieldId)
00303 {
00304   const TinyInfoOnField& info=_info_per_field[fieldId];
00305   vtkDoubleArray *arr=_arrays[_array_ids_per_field[fieldId].front()];
00306   arr->SetName(info._name.c_str());
00307   if(info._type==0)//ON_CELLS
00308     {
00309       ds->GetCellData()->AddArray(arr);
00310       return ;
00311     }
00312   if(info._type==1)//ON_NODES
00313     {
00314       ds->GetPointData()->AddArray(arr);
00315       return ;
00316     }
00317 }
00318 
00319 int ParaMEDMEM2VTK::MEDCouplingMultiFieldsFetcher::getPosGivenTimeLabel(double t)
00320 {
00321   int nbOfFields=_time_label_per_field.size();
00322   for(int i=0;i<nbOfFields;i++)
00323     if(fabs(_time_label_per_field[i]-t)<EPS_TIME)
00324       return i;
00325   //2nd chance
00326   std::vector<double>::iterator it=std::find_if(_time_label_per_field.begin(),_time_label_per_field.end(),
00327                                                 std::bind2nd(std::greater<double>(),t));
00328   if(it!=_time_label_per_field.end() && it!=_time_label_per_field.end())
00329     return std::distance(_time_label_per_field.begin(),it);
00330   //
00331   std::ostringstream oss;
00332   oss << "Unexisting time : " << t << " Not in ";
00333   std::copy(_time_label_per_field.begin(),_time_label_per_field.end(),std::ostream_iterator<double>(oss," "));
00334   oss << " !";
00335   vtkOutputWindowDisplayErrorText(oss.str().c_str());
00336   return -1;
00337 }