src/Plugins/MEDReader/IO/vtkMEDReader.cxx

   1 // Copyright (C) 2010-2015  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License, or (at your option) any later version.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay
  20
  21 #include "vtkMEDReader.h"
  22 #include "vtkGenerateVectors.h"
  23
  24 #include "vtkMultiBlockDataSet.h"
  25 #include "vtkInformation.h"
  26 #include "vtkDataSetAttributes.h"
  27 #include "vtkStringArray.h"
  28 #include "vtkMutableDirectedGraph.h"
  29 #include "vtkInformationStringKey.h"
  30 //
  31 #include "vtkUnsignedCharArray.h"
  32 #include "vtkInformationVector.h"
  33 #include "vtkSmartPointer.h"
  34 #include "vtkVariantArray.h"
  35 #include "vtkExecutive.h"
  36 #include "vtkStreamingDemandDrivenPipeline.h"
  37 #include "vtkMultiTimeStepAlgorithm.h"
  38 #include "vtkUnstructuredGrid.h"
  39 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
  40 #include "vtkQuadratureSchemeDefinition.h"
  41 #include "vtkPointData.h"
  42 #include "vtkCellData.h"
  43 #include "vtkCellType.h"
  44 #include "vtkCellArray.h"
  45 #include "vtkDoubleArray.h"
  46 #include "vtkObjectFactory.h"
  47 #ifdef MEDREADER_USE_MPI
  48 #include "vtkMultiProcessController.h"
  49 #endif
  50
  51 #include "MEDFileFieldRepresentationTree.hxx"
  52
  53 #include <string>
  54 #include <vector>
  55 #include <sstream>
  56 #include <algorithm>
  57
  58 /*!
  59  * This class stores properties in loading state mode (pvsm) when the MED file has not been read yet.
  60  * The file is not read beacause FileName has not been informed yet ! So this class stores properties of vtkMEDReader instance that
  61  * owns it and wait the vtkMEDReader::SetFileName to apply properties afterwards.
  62  */
  63 class PropertyKeeper
  64 {
  65 public:
  66   PropertyKeeper(vtkMEDReader *master):_master(master),IsGVActivated(false),GVValue(0),IsCMActivated(false),CMValue(0) { }
  67   void assignPropertiesIfNeeded();
  68   bool arePropertiesOnTreeToSetAfter() const;
  69   //
  70   void pushFieldStatusEntry(const char* name, int status);
  71   void pushGenerateVectorsValue(int value);
  72   void pushChangeModeValue(int value);
  73   void pushTimesFlagsStatusEntry(const char* name, int status);
  74 protected:
  75   // pool of pairs to assign in SetFieldsStatus if needed. The use case is the load using pvsm.
  76   std::vector< std::pair<std::string,int> > SetFieldsStatusPairs;
  77   // generate vector
  78   bool IsGVActivated;
  79   int GVValue;
  80   // change mode
  81   bool IsCMActivated;
  82   int CMValue;
  83   //
  84   std::vector< std::pair<std::string,int> > TimesFlagsStatusPairs;
  85   vtkMEDReader *_master;
  86 };
  87
  88 void PropertyKeeper::assignPropertiesIfNeeded()
  89 {
  90   if(!this->SetFieldsStatusPairs.empty())
  91     {
  92       for(std::vector< std::pair<std::string,int> >::const_iterator it=this->SetFieldsStatusPairs.begin();it!=this->SetFieldsStatusPairs.end();it++)
  93         _master->SetFieldsStatus((*it).first.c_str(),(*it).second);
  94       this->SetFieldsStatusPairs.clear();
  95     }
  96   if(!this->TimesFlagsStatusPairs.empty())
  97     {
  98       for(std::vector< std::pair<std::string,int> >::const_iterator it=this->TimesFlagsStatusPairs.begin();it!=this->TimesFlagsStatusPairs.end();it++)
  99         _master->SetTimesFlagsStatus((*it).first.c_str(),(*it).second);
 100       this->TimesFlagsStatusPairs.clear();
 101     }
 102   if(this->IsGVActivated)
 103     {
 104       _master->GenerateVectors(this->GVValue);
 105       this->IsGVActivated=false;
 106     }
 107   if(this->IsCMActivated)
 108     {
 109       _master->ChangeMode(this->CMValue);
 110       this->IsCMActivated=false;
 111     }
 112 }
 113
 114 void PropertyKeeper::pushFieldStatusEntry(const char* name, int status)
 115 {
 116   bool found(false);
 117   for(std::vector< std::pair<std::string,int> >::const_iterator it=this->SetFieldsStatusPairs.begin();it!=this->SetFieldsStatusPairs.end() && !found;it++)
 118     found=(*it).first==name;
 119   if(!found)
 120     this->SetFieldsStatusPairs.push_back(std::pair<std::string,int>(name,status));
 121 }
 122
 123 void PropertyKeeper::pushTimesFlagsStatusEntry(const char* name, int status)
 124 {
 125   bool found(false);
 126   for(std::vector< std::pair<std::string,int> >::const_iterator it=this->TimesFlagsStatusPairs.begin();it!=this->TimesFlagsStatusPairs.end() && !found;it++)
 127     found=(*it).first==name;
 128   if(!found)
 129     this->TimesFlagsStatusPairs.push_back(std::pair<std::string,int>(name,status));
 130 }
 131
 132 void PropertyKeeper::pushGenerateVectorsValue(int value)
 133 {
 134   this->IsGVActivated=true;
 135   this->GVValue=value;
 136 }
 137
 138 void PropertyKeeper::pushChangeModeValue(int value)
 139 {
 140   this->IsCMActivated=true;
 141   this->CMValue=value;
 142 }
 143
 144 bool PropertyKeeper::arePropertiesOnTreeToSetAfter() const
 145 {
 146   return !SetFieldsStatusPairs.empty();
 147 }
 148
 149 class vtkMEDReader::vtkMEDReaderInternal
 150 {
 151
 152 public:
 153   vtkMEDReaderInternal(vtkMEDReader *master):TK(0),IsMEDOrSauv(true),IsStdOrMode(false),GenerateVect(false),SIL(0),LastLev0(-1),FirstCall0(2),PK(master),MyMTime(0)
 154   {
 155   }
 156
 157   bool PluginStart0()
 158   {
 159     if(FirstCall0==0)
 160       return false;
 161     FirstCall0--;
 162     return true;
 163   }
 164
 165   ~vtkMEDReaderInternal()
 166   {
 167     if(this->SIL)
 168       this->SIL->Delete();
 169   }
 170 public:
 171   MEDFileFieldRepresentationTree Tree;
 172   TimeKeeper TK;
 173   std::string FileName;
 174   //when true the file is MED file. when false it is a Sauv file
 175   bool IsMEDOrSauv;
 176   //when false -> std, true -> mode. By default std (false).
 177   bool IsStdOrMode;
 178   //when false -> do nothing. When true cut off or extend to nbOfCompo=3 vector arrays.
 179   bool GenerateVect;
 180   std::string DftMeshName;
 181   // Store the vtkMutableDirectedGraph that represents links between family, groups and cell types
 182   vtkMutableDirectedGraph* SIL;
 183   // store the lev0 id in Tree corresponding to the TIME_STEPS in the pipeline.
 184   int LastLev0;
 185   // The property keeper is usable only in pvsm mode.
 186   PropertyKeeper PK;
 187   int MyMTime;
 188   std::set<std::string> _wonderful_set;// this set is used by SetFieldsStatus method to detect the fact that SetFieldsStatus has been called for all items ! Great !
 189 private:
 190   unsigned char FirstCall0;
 191 };
 192
 193 vtkStandardNewMacro(vtkMEDReader);
 194
 195 vtkMEDReader::vtkMEDReader():Internal(new vtkMEDReaderInternal(this))
 196 {
 197   this->SetNumberOfInputPorts(0);
 198   this->SetNumberOfOutputPorts(1);
 199 }
 200
 201 vtkMEDReader::~vtkMEDReader()
 202 {
 203   delete this->Internal;
 204 }
 205
 206 void vtkMEDReader::Reload(int a)
 207 {
 208   if(a==0)
 209     return;
 210   std::cerr << "vtkMEDReader::Reload" << a << std::endl;
 211   std::string fName((const char *)this->GetFileName());
 212   delete this->Internal;
 213   this->Internal=new vtkMEDReaderInternal(this);
 214   this->SetFileName(fName.c_str());
 215 }
 216
 217 int vtkMEDReader::GetServerModifTime()
 218 {
 219   return this->Internal->MyMTime;
 220 }
 221
 222 void vtkMEDReader::GenerateVectors(int val)
 223 {
 224   if(this->Internal->FileName.empty())
 225     {//pvsm mode
 226       this->Internal->PK.pushGenerateVectorsValue(val);
 227       return ;
 228     }
 229   //not pvsm mode (general case)
 230   bool val2((bool)val);
 231   if(val2!=this->Internal->GenerateVect)
 232     {
 233       this->Internal->GenerateVect=val2;
 234       this->Modified();
 235     }
 236 }
 237
 238 void vtkMEDReader::ChangeMode(int newMode)
 239 {
 240   if(this->Internal->FileName.empty())
 241     {//pvsm mode
 242       this->Internal->PK.pushChangeModeValue(newMode);
 243       return ;
 244     }
 245   //not pvsm mode (general case)
 246   this->Internal->IsStdOrMode=newMode!=0;
 247   //std::cerr << "vtkMEDReader::ChangeMode : " << this->Internal->IsStdOrMode << std::endl;
 248   this->Modified();
 249 }
 250
 251 const char *vtkMEDReader::GetSeparator()
 252 {
 253   return MEDFileFieldRepresentationLeavesArrays::ZE_SEP;
 254 }
 255
 256 void vtkMEDReader::SetFileName(const char *fname)
 257 {
 258   try
 259     {
 260       this->Internal->FileName=fname;
 261       std::size_t pos(this->Internal->FileName.find_last_of('.'));
 262       if(pos!=std::string::npos)
 263         {
 264           std::string ext(this->Internal->FileName.substr(pos));
 265           if(ext.find("sauv")!=std::string::npos)
 266             this->Internal->IsMEDOrSauv=false;
 267         }
 268       if(this->Internal->Tree.getNumberOfLeavesArrays()==0)
 269         {
 270           int iPart(-1),nbOfParts(-1);
 271 #ifdef MEDREADER_USE_MPI
 272           vtkMultiProcessController *vmpc(vtkMultiProcessController::GetGlobalController());
 273           if(vmpc)
 274             {
 275               iPart=vmpc->GetLocalProcessId();
 276               nbOfParts=vmpc->GetNumberOfProcesses();
 277             }
 278 #endif
 279           this->Internal->Tree.loadMainStructureOfFile(this->Internal->FileName.c_str(),this->Internal->IsMEDOrSauv,iPart,nbOfParts);
 280           if(!this->Internal->PK.arePropertiesOnTreeToSetAfter())
 281             this->Internal->Tree.activateTheFirst();//This line manually initialize the status of server (this) with the remote client.
 282           this->Internal->TK.setMaxNumberOfTimeSteps(this->Internal->Tree.getMaxNumberOfTimeSteps());
 283         }
 284       this->Modified();
 285       this->Internal->PK.assignPropertiesIfNeeded();
 286     }
 287   catch(INTERP_KERNEL::Exception& e)
 288     {
 289       delete this->Internal;
 290       this->Internal=0;
 291       std::ostringstream oss;
 292       oss << "Exception has been thrown in vtkMEDReader::SetFileName : " << e.what() << std::endl;
 293       if(this->HasObserver("ErrorEvent") )
 294         this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
 295       else
 296         vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
 297       vtkObject::BreakOnError();
 298     }
 299 }
 300
 301 char *vtkMEDReader::GetFileName()
 302 {
 303   return const_cast<char *>(this->Internal->FileName.c_str());
 304 }
 305
 306 int vtkMEDReader::RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 307 {
 308   //std::cerr << "########################################## RequestInformation ##########################################" << std::endl;
 309   if(!this->Internal)
 310     return 0;
 311   try
 312     {
 313       vtkInformation *outInfo(outputVector->GetInformationObject(0));
 314       outInfo->Set(vtkDataObject::DATA_TYPE_NAME(),"vtkMultiBlockDataSet");
 315       this->UpdateSIL(outInfo);
 316       //
 317       bool dummy(false);
 318       this->PublishTimeStepsIfNeeded(outInfo,dummy);
 319     }
 320   catch(INTERP_KERNEL::Exception& e)
 321     {
 322       std::ostringstream oss;
 323       oss << "Exception has been thrown in vtkMEDReader::RequestInformation : " << e.what() << std::endl;
 324       if(this->HasObserver("ErrorEvent") )
 325         this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
 326       else
 327         vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
 328       vtkObject::BreakOnError();
 329       return 0;
 330     }
 331   return 1;
 332 }
 333
 334 int vtkMEDReader::RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 335 {
 336   //std::cerr << "########################################## RequestData        ##########################################";
 337   if(!this->Internal)
 338     return 0;
 339   try
 340     {
 341       vtkInformation *outInfo(outputVector->GetInformationObject(0));
 342       vtkMultiBlockDataSet *output(vtkMultiBlockDataSet::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
 343       bool isUpdated(false);
 344       double reqTS(0.);
 345       if(outInfo->Has(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP()))
 346         reqTS=outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP());
 347       //std::cerr << reqTS << std::endl;
 348       this->FillMultiBlockDataSetInstance(output,reqTS);
 349       output->GetInformation()->Set(vtkDataObject::DATA_TIME_STEP(),reqTS);
 350       this->UpdateSIL(outInfo);
 351       //this->UpdateProgress((float) progress/((float) maxprogress-1));
 352     }
 353   catch(INTERP_KERNEL::Exception& e)
 354     {
 355       std::cerr << "Exception has been thrown in vtkMEDReader::RequestData : " << e.what() << std::endl;
 356       return 0;
 357     }
 358   return 1;
 359 }
 360
 361 void vtkMEDReader::SetFieldsStatus(const char* name, int status)
 362 {
 363   //this->Internal->_wonderful_set.insert(name);
 364   if(this->Internal->FileName.empty())
 365     {//pvsm mode
 366       this->Internal->PK.pushFieldStatusEntry(name,status);
 367       return ;
 368     }
 369   this->Internal->_wonderful_set.insert(name);
 370   //not pvsm mode (general case)
 371   try
 372     {
 373       this->Internal->Tree.changeStatusOfAndUpdateToHaveCoherentVTKDataSet(this->Internal->Tree.getIdHavingZeName(name),status);
 374       if(this->Internal->_wonderful_set.size()==GetNumberOfFieldsTreeArrays())
 375         {
 376           if(!this->Internal->PluginStart0())
 377             {
 378               this->Modified();
 379             }
 380           this->Internal->MyMTime++;
 381           this->Internal->_wonderful_set.clear();
 382         }
 383     }
 384   catch(INTERP_KERNEL::Exception& e)
 385     {
 386       if(!this->Internal->FileName.empty())
 387         {
 388           std::cerr << "vtkMEDReader::SetFieldsStatus error : " << e.what() << " *** WITH STATUS=" << status << endl;
 389           return ;
 390         }
 391     }
 392 }
 393
 394 int vtkMEDReader::GetNumberOfFieldsTreeArrays()
 395 {
 396   if(!this->Internal)
 397     return 0;
 398   return this->Internal->Tree.getNumberOfLeavesArrays();
 399   //std::cerr << "vtkMEDReader::GetNumberOfFieldsTreeArrays called ! " << ret << std::endl;
 400 }
 401
 402 const char *vtkMEDReader::GetFieldsTreeArrayName(int index)
 403 {
 404   return this->Internal->Tree.getNameOfC(index);
 405   //std::cerr << "vtkMEDReader::GetFieldsTreeArrayName(" << index << ") called ! " << ret << std::endl;
 406 }
 407
 408 int vtkMEDReader::GetFieldsTreeArrayStatus(const char *name)
 409 {
 410   int zeId(this->Internal->Tree.getIdHavingZeName(name));
 411   int ret(this->Internal->Tree.getStatusOf(zeId));
 412   return ret;
 413 }
 414
 415 void vtkMEDReader::SetTimesFlagsStatus(const char *name, int status)
 416 {
 417   if(this->Internal->FileName.empty())
 418     {//pvsm mode
 419       this->Internal->PK.pushTimesFlagsStatusEntry(name,status);
 420       return ;
 421     }
 422   //not pvsm mode (general case)
 423   int pos(0);
 424   std::istringstream iss(name); iss >> pos;
 425   this->Internal->TK.getTimesFlagArray()[pos].first=(bool)status;
 426   if(pos==this->Internal->TK.getTimesFlagArray().size()-1)
 427     if(!this->Internal->PluginStart0())
 428       {
 429         this->Modified();
 430         //this->Internal->TK.printSelf(std::cerr);
 431       }
 432 }
 433
 434 int vtkMEDReader::GetNumberOfTimesFlagsArrays()
 435 {
 436   if(!this->Internal)
 437     return 0;
 438   return (int)this->Internal->TK.getTimesFlagArray().size();
 439 }
 440
 441 const char *vtkMEDReader::GetTimesFlagsArrayName(int index)
 442 {
 443   return this->Internal->TK.getTimesFlagArray()[index].second.c_str();
 444 }
 445
 446 int vtkMEDReader::GetTimesFlagsArrayStatus(const char *name)
 447 {
 448   int pos(0);
 449   std::istringstream iss(name); iss >> pos;
 450   return (int)this->Internal->TK.getTimesFlagArray()[pos].first;
 451 }
 452
 453 void vtkMEDReader::UpdateSIL(vtkInformation *info)
 454 {
 455   if(!this->Internal)
 456       return ;
 457   vtkMutableDirectedGraph *sil(vtkMutableDirectedGraph::New());
 458   std::string meshName(this->BuildSIL(sil));
 459   if(meshName!=this->Internal->DftMeshName)
 460     {
 461       if(this->Internal->SIL)
 462         this->Internal->SIL->Delete();
 463       this->Internal->SIL=sil;
 464       this->Internal->DftMeshName=meshName;
 465       info->Set(vtkDataObject::SIL(),this->Internal->SIL);
 466       //request->AppendUnique(vtkExecutive::KEYS_TO_COPY(),vtkDataObject::SIL());
 467     }
 468   else
 469     {
 470       sil->Delete();
 471     }
 472 }
 473
 474 /*!
 475  * The returned string is the name of the mesh activated which groups and families are in \a sil.
 476  */
 477 std::string vtkMEDReader::BuildSIL(vtkMutableDirectedGraph* sil)
 478 {
 479   sil->Initialize();
 480   vtkSmartPointer<vtkVariantArray> childEdge(vtkSmartPointer<vtkVariantArray>::New());
 481   childEdge->InsertNextValue(0);
 482   vtkSmartPointer<vtkVariantArray> crossEdge(vtkSmartPointer<vtkVariantArray>::New());
 483   crossEdge->InsertNextValue(1);
 484   // CrossEdge is an edge linking hierarchies.
 485   vtkUnsignedCharArray* crossEdgesArray=vtkUnsignedCharArray::New();
 486   crossEdgesArray->SetName("CrossEdges");
 487   sil->GetEdgeData()->AddArray(crossEdgesArray);
 488   crossEdgesArray->Delete();
 489   std::vector<std::string> names;
 490   // Now build the hierarchy.
 491   vtkIdType rootId=sil->AddVertex();
 492   names.push_back("SIL");
 493   // Add global fields root
 494   vtkIdType fieldsRoot(sil->AddChild(rootId,childEdge));
 495   names.push_back("FieldsStatusTree");
 496   this->Internal->Tree.feedSIL(sil,fieldsRoot,childEdge,names);
 497   vtkIdType meshesFamsGrpsRoot(sil->AddChild(rootId,childEdge));
 498   names.push_back("MeshesFamsGrps");
 499   std::string dftMeshName(this->Internal->Tree.feedSILForFamsAndGrps(sil,meshesFamsGrpsRoot,childEdge,names));
 500   // This array is used to assign names to nodes.
 501   vtkStringArray *namesArray(vtkStringArray::New());
 502   namesArray->SetName("Names");
 503   namesArray->SetNumberOfTuples(sil->GetNumberOfVertices());
 504   sil->GetVertexData()->AddArray(namesArray);
 505   namesArray->Delete();
 506   std::vector<std::string>::const_iterator iter;
 507   vtkIdType cc;
 508   for(cc=0, iter=names.begin(); iter!=names.end(); ++iter, ++cc)
 509     namesArray->SetValue(cc,(*iter).c_str());
 510   return dftMeshName;
 511 }
 512
 513 double vtkMEDReader::PublishTimeStepsIfNeeded(vtkInformation *outInfo, bool& isUpdated)
 514 {
 515   int lev0(-1);
 516   std::vector<double> tsteps;
 517   if(!this->Internal->IsStdOrMode)
 518     tsteps=this->Internal->Tree.getTimeSteps(lev0,this->Internal->TK);
 519   else
 520     { tsteps.resize(1); tsteps[0]=0.; }
 521   isUpdated=false;
 522   if(lev0!=this->Internal->LastLev0)
 523     {
 524       isUpdated=true;
 525       double timeRange[2];
 526       timeRange[0]=tsteps.front();
 527       timeRange[1]=tsteps.back();
 528       outInfo->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(),&tsteps[0],tsteps.size());
 529       outInfo->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(),timeRange,2);
 530       this->Internal->LastLev0=lev0;
 531     }
 532   return tsteps.front();
 533 }
 534
 535 void vtkMEDReader::FillMultiBlockDataSetInstance(vtkMultiBlockDataSet *output, double reqTS)
 536 {
 537   std::string meshName;
 538   vtkDataSet *ret(this->Internal->Tree.buildVTKInstance(this->Internal->IsStdOrMode,reqTS,meshName,this->Internal->TK));
 539   if(this->Internal->GenerateVect)
 540     {
 541       vtkGenerateVectors::Operate(ret->GetPointData());
 542       vtkGenerateVectors::Operate(ret->GetCellData());
 543       vtkGenerateVectors::Operate(ret->GetFieldData());
 544     }
 545   output->SetBlock(0,ret);
 546   ret->Delete();
 547 }
 548
 549 void vtkMEDReader::PrintSelf(ostream& os, vtkIndent indent)
 550 {
 551   this->Superclass::PrintSelf(os, indent);
 552 }