1 // Copyright (C) 2010-2020 CEA/DEN, EDF R&D
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.
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.
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
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay
21 #include "vtkMEDReader.h"
22 #include "vtkGenerateVectors.h"
23 #include "MEDUtilities.hxx"
25 #include "vtkMultiBlockDataSet.h"
26 #include "vtkInformation.h"
27 #include "vtkDataSetAttributes.h"
28 #include "vtkStringArray.h"
29 #include "vtkMutableDirectedGraph.h"
30 #include "vtkInformationStringKey.h"
32 #include "vtkUnsignedCharArray.h"
33 #include "vtkInformationVector.h"
34 #include "vtkSmartPointer.h"
35 #include "vtkVariantArray.h"
36 #include "vtkExecutive.h"
37 #include "vtkStreamingDemandDrivenPipeline.h"
38 #include "vtkMultiTimeStepAlgorithm.h"
39 #include "vtkUnstructuredGrid.h"
40 #include "vtkInformationQuadratureSchemeDefinitionVectorKey.h"
41 #include "vtkInformationDoubleVectorKey.h"
42 #include "vtkQuadratureSchemeDefinition.h"
43 #include "vtkPointData.h"
44 #include "vtkCellData.h"
45 #include "vtkCellType.h"
46 #include "vtkCellArray.h"
47 #include "vtkDoubleArray.h"
48 #include "vtkObjectFactory.h"
49 #include "vtkInformationDataObjectMetaDataKey.h"
51 #ifdef MEDREADER_USE_MPI
52 #include "vtkMultiProcessController.h"
53 #include "vtkPUnstructuredGridGhostCellsGenerator.h"
56 #include "MEDFileFieldRepresentationTree.hxx"
65 * This class stores properties in loading state mode (pvsm) when the MED file has not been read yet.
66 * The file is not read beacause FileName has not been informed yet ! So this class stores properties of vtkMEDReader instance that
67 * owns it and wait the vtkMEDReader::SetFileName to apply properties afterwards.
72 PropertyKeeper(vtkMEDReader *master):IsGVActivated(false),GVValue(0),IsCMActivated(false),CMValue(0),IsGhostActivated(false),GCGCP(1),_master(master) { }
73 void assignPropertiesIfNeeded();
74 bool arePropertiesOnTreeToSetAfter() const;
76 void pushFieldStatusEntry(const char* name, int status);
77 void pushGenerateVectorsValue(int value);
78 void pushChangeModeValue(int value);
79 void pushTimesFlagsStatusEntry(const char* name, int status);
80 void pushGhost(int value);
82 // pool of pairs to assign in SetFieldsStatus if needed. The use case is the load using pvsm.
83 std::vector< std::pair<std::string,int> > SetFieldsStatusPairs;
91 bool IsGhostActivated;
94 std::vector< std::pair<std::string,int> > TimesFlagsStatusPairs;
95 vtkMEDReader *_master;
98 void PropertyKeeper::assignPropertiesIfNeeded()
100 if(!this->SetFieldsStatusPairs.empty())
102 for(std::vector< std::pair<std::string,int> >::const_iterator it=this->SetFieldsStatusPairs.begin();it!=this->SetFieldsStatusPairs.end();it++)
103 _master->SetFieldsStatus((*it).first.c_str(),(*it).second);
104 this->SetFieldsStatusPairs.clear();
106 if(!this->TimesFlagsStatusPairs.empty())
108 for(std::vector< std::pair<std::string,int> >::const_iterator it=this->TimesFlagsStatusPairs.begin();it!=this->TimesFlagsStatusPairs.end();it++)
109 _master->SetTimesFlagsStatus((*it).first.c_str(),(*it).second);
110 this->TimesFlagsStatusPairs.clear();
112 if(this->IsGVActivated)
114 _master->GenerateVectors(this->GVValue);
115 this->IsGVActivated=false;
117 if(this->IsCMActivated)
119 _master->ChangeMode(this->CMValue);
120 this->IsCMActivated=false;
122 if(this->IsGhostActivated)
124 _master->GhostCellGeneratorCallForPara(this->GCGCP);
125 this->IsGhostActivated=false;
129 void PropertyKeeper::pushFieldStatusEntry(const char* name, int status)
132 for(std::vector< std::pair<std::string,int> >::const_iterator it=this->SetFieldsStatusPairs.begin();it!=this->SetFieldsStatusPairs.end() && !found;it++)
133 found=(*it).first==name;
135 this->SetFieldsStatusPairs.push_back(std::pair<std::string,int>(name,status));
138 void PropertyKeeper::pushTimesFlagsStatusEntry(const char* name, int status)
141 for(std::vector< std::pair<std::string,int> >::const_iterator it=this->TimesFlagsStatusPairs.begin();it!=this->TimesFlagsStatusPairs.end() && !found;it++)
142 found=(*it).first==name;
144 this->TimesFlagsStatusPairs.push_back(std::pair<std::string,int>(name,status));
147 void PropertyKeeper::pushGenerateVectorsValue(int value)
149 this->IsGVActivated=true;
153 void PropertyKeeper::pushChangeModeValue(int value)
155 this->IsCMActivated=true;
159 void PropertyKeeper::pushGhost(int value)
161 this->IsGhostActivated=true;
165 bool PropertyKeeper::arePropertiesOnTreeToSetAfter() const
167 return !SetFieldsStatusPairs.empty();
170 class vtkMEDReader::vtkMEDReaderInternal
174 vtkMEDReaderInternal(vtkMEDReader *master):TK(0),IsMEDOrSauv(true),IsStdOrMode(false),GenerateVect(false),SIL(0),LastLev0(-1),PK(master),MyMTime(0),GCGCP(true),FirstCall0(2)
180 return false; // TODO Useless and buggy
187 ~vtkMEDReaderInternal()
193 MEDFileFieldRepresentationTree Tree;
195 std::string FileName;
196 //when true the file is MED file. when false it is a Sauv file
198 //when false -> std, true -> mode. By default std (false).
200 //when false -> do nothing. When true cut off or extend to nbOfCompo=3 vector arrays.
202 std::string DftMeshName;
203 // Store the vtkMutableDirectedGraph that represents links between family, groups and cell types
204 vtkMutableDirectedGraph* SIL;
205 // store the lev0 id in Tree corresponding to the TIME_STEPS in the pipeline.
207 // The property keeper is usable only in pvsm mode.
210 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 Items are not sorted ! Why ?
211 std::map<std::string,bool> _wonderful_ref;// this map stores the state before a SetFieldsStatus status.
215 unsigned char FirstCall0;
218 vtkStandardNewMacro(vtkMEDReader)
220 // vtkInformationKeyMacro(vtkMEDReader, META_DATA, DataObjectMetaData); // Here we need to customize vtkMEDReader::META_DATA method
221 // start of overload of vtkInformationKeyMacro
222 static vtkInformationDataObjectMetaDataKey *vtkMEDReader_META_DATA=new vtkInformationDataObjectMetaDataKey("META_DATA","vtkMEDReader");
224 vtkInformationDataObjectMetaDataKey *vtkMEDReader::META_DATA()
226 static const char ZE_KEY[]="vtkMEDReader::META_DATA";
227 vtkInformationDataObjectMetaDataKey *ret(vtkMEDReader_META_DATA);
228 MEDCoupling::GlobalDict *gd(MEDCoupling::GlobalDict::GetInstance());
229 if(!gd->hasKey(ZE_KEY))
230 {// here META_DATA is put on global var to be exchanged with other filters without dependancy of MEDReader. Please do not change ZE_KEY !
231 std::ostringstream oss; oss << ret;
232 gd->setKeyValue(ZE_KEY,oss.str());
237 static vtkInformationGaussDoubleVectorKey *vtkMEDReader_GAUSS_DATA=new vtkInformationGaussDoubleVectorKey("GAUSS_DATA","vtkMEDReader");
239 vtkInformationGaussDoubleVectorKey *vtkMEDReader::GAUSS_DATA()
241 static const char ZE_KEY[]="vtkMEDReader::GAUSS_DATA";
242 vtkInformationGaussDoubleVectorKey *ret(vtkMEDReader_GAUSS_DATA);
243 MEDCoupling::GlobalDict *gd(MEDCoupling::GlobalDict::GetInstance());
244 if(!gd->hasKey(ZE_KEY))
245 {// here META_DATA is put on global var to be exchanged with other filters without dependancy of MEDReader. Please do not change ZE_KEY !
246 vtkInformationDoubleVectorKey *ret2(ret);
247 std::ostringstream oss; oss << ret2;
248 gd->setKeyValue(ZE_KEY,oss.str());
252 // end of overload of vtkInformationKeyMacro
254 vtkMEDReader::vtkMEDReader():Internal(new vtkMEDReaderInternal(this))
256 this->SetNumberOfInputPorts(0);
257 this->SetNumberOfOutputPorts(1);
260 vtkMEDReader::~vtkMEDReader()
262 delete this->Internal;
266 void vtkMEDReader::Reload()
268 std::string fName((const char *)this->GetFileName());
269 delete this->Internal;
270 this->Internal=new vtkMEDReaderInternal(this);
271 this->SetFileName(fName.c_str());
274 int vtkMEDReader::GetServerModifTime()
276 if( !this->Internal )
278 return this->Internal->MyMTime;
281 void vtkMEDReader::GenerateVectors(int val)
283 if ( !this->Internal )
286 if(this->Internal->FileName.empty())
288 this->Internal->PK.pushGenerateVectorsValue(val);
291 //not pvsm mode (general case)
292 bool val2((bool)val);
293 if(val2!=this->Internal->GenerateVect)
295 this->Internal->GenerateVect=val2;
300 void vtkMEDReader::ChangeMode(int newMode)
302 if ( !this->Internal )
305 if(this->Internal->FileName.empty())
307 this->Internal->PK.pushChangeModeValue(newMode);
310 //not pvsm mode (general case)
311 this->Internal->IsStdOrMode=newMode!=0;
315 void vtkMEDReader::GhostCellGeneratorCallForPara(int gcgcp)
317 if ( !this->Internal )
320 if(this->Internal->FileName.empty())
322 this->Internal->PK.pushGhost(gcgcp);
325 bool newVal(gcgcp!=0);
326 if(newVal!=this->Internal->GCGCP)
328 this->Internal->GCGCP=newVal;
333 const char *vtkMEDReader::GetSeparator()
335 return MEDFileFieldRepresentationLeavesArrays::ZE_SEP;
338 void vtkMEDReader::SetFileName(const char *fname)
344 this->Internal->FileName=fname;
345 std::size_t pos(this->Internal->FileName.find_last_of('.'));
346 if(pos!=std::string::npos)
348 std::string ext(this->Internal->FileName.substr(pos));
349 if(ext.find("sauv")!=std::string::npos)
350 this->Internal->IsMEDOrSauv=false;
352 if(this->Internal->Tree.getNumberOfLeavesArrays()==0)
354 int iPart(-1),nbOfParts(-1);
355 #ifdef MEDREADER_USE_MPI
356 vtkMultiProcessController *vmpc(vtkMultiProcessController::GetGlobalController());
359 iPart=vmpc->GetLocalProcessId();
360 nbOfParts=vmpc->GetNumberOfProcesses();
363 this->Internal->Tree.loadMainStructureOfFile(this->Internal->FileName.c_str(),this->Internal->IsMEDOrSauv,iPart,nbOfParts);
364 if(!this->Internal->PK.arePropertiesOnTreeToSetAfter())
365 this->Internal->Tree.activateTheFirst();//This line manually initialize the status of server (this) with the remote client.
366 this->Internal->TK.setMaxNumberOfTimeSteps(this->Internal->Tree.getMaxNumberOfTimeSteps());
369 this->Internal->PK.assignPropertiesIfNeeded();
371 catch(INTERP_KERNEL::Exception& e)
373 delete this->Internal;
375 std::ostringstream oss;
376 oss << "Exception has been thrown in vtkMEDReader::SetFileName : " << e.what() << std::endl;
377 if(this->HasObserver("ErrorEvent") )
378 this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
380 vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
381 vtkObject::BreakOnError();
385 char *vtkMEDReader::GetFileName()
389 return const_cast<char *>(this->Internal->FileName.c_str());
392 int vtkMEDReader::RequestInformation(vtkInformation *request, vtkInformationVector ** /*inputVector*/, vtkInformationVector *outputVector)
394 // std::cout << "########################################## vtkMEDReader::RequestInformation ##########################################" << std::endl;
399 // request->Print(cout);
400 vtkInformation *outInfo(outputVector->GetInformationObject(0));
401 outInfo->Set(vtkDataObject::DATA_TYPE_NAME(),"vtkMultiBlockDataSet");
402 this->UpdateSIL(request, outInfo);
404 // Set the meta data graph as a meta data key in the information
405 // That's all that is needed to transfer it along the pipeline
406 outInfo->Set(vtkMEDReader::META_DATA(),this->Internal->SIL);
409 this->PublishTimeStepsIfNeeded(outInfo,dummy);
411 catch(INTERP_KERNEL::Exception& e)
413 std::ostringstream oss;
414 oss << "Exception has been thrown in vtkMEDReader::RequestInformation : " << e.what() << std::endl;
415 if(this->HasObserver("ErrorEvent") )
416 this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
418 vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
419 vtkObject::BreakOnError();
425 int vtkMEDReader::RequestData(vtkInformation *request, vtkInformationVector ** /*inputVector*/, vtkInformationVector *outputVector)
427 // std::cout << "########################################## vtkMEDReader::RequestData ##########################################" << std::endl;
432 // request->Print(cout);
433 vtkInformation *outInfo(outputVector->GetInformationObject(0));
434 vtkMultiBlockDataSet *output(vtkMultiBlockDataSet::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
435 //bool isUpdated(false); // todo: unused
437 if(outInfo->Has(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP()))
438 reqTS=outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP());
440 #ifndef MEDREADER_USE_MPI
441 this->FillMultiBlockDataSetInstance(output,reqTS,&ti);
443 if(this->Internal->GCGCP)
445 vtkSmartPointer<vtkPUnstructuredGridGhostCellsGenerator> gcg(vtkSmartPointer<vtkPUnstructuredGridGhostCellsGenerator>::New());
447 vtkDataSet *ret(RetrieveDataSetAtTime(reqTS,&ti));
448 gcg->SetInputData(ret);
451 gcg->SetUseGlobalPointIds(true);
452 gcg->SetBuildIfRequired(false);
454 output->SetBlock(0,gcg->GetOutput());
457 this->FillMultiBlockDataSetInstance(output,reqTS,&ti);
461 const std::vector<double>& data(ti.getData());
462 outInfo->Set(vtkMEDReader::GAUSS_DATA(),&data[0],(int)data.size());
463 request->Append(vtkExecutive::KEYS_TO_COPY(),vtkMEDReader::GAUSS_DATA());// Thank you to SciberQuest and DIPOLE_CENTER ! Don't understand why ! In RequestInformation it does not work !
465 output->GetInformation()->Set(vtkDataObject::DATA_TIME_STEP(),reqTS);
466 // Is it really needed ? TODO
467 this->UpdateSIL(request, outInfo);
469 catch(INTERP_KERNEL::Exception& e)
471 std::cerr << "Exception has been thrown in vtkMEDReader::RequestData : " << e.what() << std::endl;
477 void vtkMEDReader::SetFieldsStatus(const char* name, int status)
479 if( !this->Internal )
482 //this->Internal->_wonderful_set.insert(name);
483 if(this->Internal->FileName.empty())
485 this->Internal->PK.pushFieldStatusEntry(name,status);
488 if(this->Internal->_wonderful_set.empty())
489 this->Internal->_wonderful_ref=this->Internal->Tree.dumpState();// start of SetFieldsStatus serie -> store ref to compare at the end of the SetFieldsStatus serie.
490 this->Internal->_wonderful_set.insert(name);
491 //not pvsm mode (general case)
494 this->Internal->Tree.changeStatusOfAndUpdateToHaveCoherentVTKDataSet(this->Internal->Tree.getIdHavingZeName(name),status);
495 if((int)this->Internal->_wonderful_set.size()==GetNumberOfFieldsTreeArrays())
497 if(this->Internal->_wonderful_ref!=this->Internal->Tree.dumpState())
499 if(!this->Internal->PluginStart0())
503 this->Internal->MyMTime++;
505 this->Internal->_wonderful_set.clear();
508 catch(INTERP_KERNEL::Exception& e)
510 if(!this->Internal->FileName.empty())
512 std::cerr << "vtkMEDReader::SetFieldsStatus error : " << e.what() << " *** WITH STATUS=" << status << endl;
518 int vtkMEDReader::GetNumberOfFieldsTreeArrays()
522 return this->Internal->Tree.getNumberOfLeavesArrays();
525 const char *vtkMEDReader::GetFieldsTreeArrayName(int index)
529 return this->Internal->Tree.getNameOfC(index);
532 int vtkMEDReader::GetFieldsTreeArrayStatus(const char *name)
537 int zeId(this->Internal->Tree.getIdHavingZeName(name));
538 int ret(this->Internal->Tree.getStatusOf(zeId));
542 void vtkMEDReader::SetTimesFlagsStatus(const char *name, int status)
547 if(this->Internal->FileName.empty())
549 this->Internal->PK.pushTimesFlagsStatusEntry(name,status);
552 //not pvsm mode (general case)
554 std::istringstream iss(name); iss >> pos;
555 this->Internal->TK.getTimesFlagArray()[pos].first=(bool)status;
556 if(pos==(int)this->Internal->TK.getTimesFlagArray().size()-1)
557 if(!this->Internal->PluginStart0())
560 //this->Internal->TK.printSelf(std::cerr);
564 int vtkMEDReader::GetNumberOfTimesFlagsArrays()
568 return (int)this->Internal->TK.getTimesFlagArray().size();
571 const char *vtkMEDReader::GetTimesFlagsArrayName(int index)
573 return this->Internal->TK.getTimesFlagArray()[index].second.c_str();
576 int vtkMEDReader::GetTimesFlagsArrayStatus(const char *name)
581 std::istringstream iss(name); iss >> pos;
582 return (int)this->Internal->TK.getTimesFlagArray()[pos].first;
585 void vtkMEDReader::UpdateSIL(vtkInformation* /*request*/, vtkInformation * /*info*/)
589 std::string meshName(this->Internal->Tree.getActiveMeshName());
590 if(!this->Internal->SIL || meshName!=this->Internal->DftMeshName)
592 vtkMutableDirectedGraph *sil(vtkMutableDirectedGraph::New());
594 if(this->Internal->SIL)
595 this->Internal->SIL->Delete();
596 this->Internal->SIL=sil;
597 this->Internal->DftMeshName=meshName;
602 * The returned string is the name of the mesh activated which groups and families are in \a sil.
604 std::string vtkMEDReader::BuildSIL(vtkMutableDirectedGraph* sil)
607 return std::string();
609 vtkSmartPointer<vtkVariantArray> childEdge(vtkSmartPointer<vtkVariantArray>::New());
610 childEdge->InsertNextValue(0);
611 vtkSmartPointer<vtkVariantArray> crossEdge(vtkSmartPointer<vtkVariantArray>::New());
612 crossEdge->InsertNextValue(1);
613 // CrossEdge is an edge linking hierarchies.
614 vtkUnsignedCharArray* crossEdgesArray=vtkUnsignedCharArray::New();
615 crossEdgesArray->SetName("CrossEdges");
616 sil->GetEdgeData()->AddArray(crossEdgesArray);
617 crossEdgesArray->Delete();
618 std::vector<std::string> names;
619 // Now build the hierarchy.
620 vtkIdType rootId=sil->AddVertex();
621 names.push_back("SIL");
622 // Add global fields root
623 vtkIdType fieldsRoot(sil->AddChild(rootId,childEdge));
624 names.push_back("FieldsStatusTree");
625 this->Internal->Tree.feedSIL(sil,fieldsRoot,childEdge,names);
626 vtkIdType meshesFamsGrpsRoot(sil->AddChild(rootId,childEdge));
627 names.push_back("MeshesFamsGrps");
628 std::string dftMeshName(this->Internal->Tree.feedSILForFamsAndGrps(sil,meshesFamsGrpsRoot,childEdge,names));
629 // This array is used to assign names to nodes.
630 vtkStringArray *namesArray(vtkStringArray::New());
631 namesArray->SetName("Names");
632 namesArray->SetNumberOfTuples(sil->GetNumberOfVertices());
633 sil->GetVertexData()->AddArray(namesArray);
634 namesArray->Delete();
635 std::vector<std::string>::const_iterator iter;
637 for(cc=0, iter=names.begin(); iter!=names.end(); ++iter, ++cc)
638 namesArray->SetValue(cc,(*iter).c_str());
642 double vtkMEDReader::PublishTimeStepsIfNeeded(vtkInformation *outInfo, bool& isUpdated)
648 std::vector<double> tsteps;
649 if(!this->Internal->IsStdOrMode)
650 tsteps=this->Internal->Tree.getTimeSteps(lev0,this->Internal->TK);
652 { tsteps.resize(1); tsteps[0]=0.; }
654 if(lev0!=this->Internal->LastLev0)
658 timeRange[0]=tsteps.front();
659 timeRange[1]=tsteps.back();
660 outInfo->Set(vtkStreamingDemandDrivenPipeline::TIME_STEPS(),&tsteps[0],(int)tsteps.size());
661 outInfo->Set(vtkStreamingDemandDrivenPipeline::TIME_RANGE(),timeRange,2);
662 this->Internal->LastLev0=lev0;
664 return tsteps.front();
667 void vtkMEDReader::FillMultiBlockDataSetInstance(vtkMultiBlockDataSet *output, double reqTS, ExportedTinyInfo *internalInfo)
669 if( !this->Internal )
671 vtkDataSet *ret(RetrieveDataSetAtTime(reqTS,internalInfo));
672 output->SetBlock(0,ret);
676 vtkDataSet *vtkMEDReader::RetrieveDataSetAtTime(double reqTS, ExportedTinyInfo *internalInfo)
678 if( !this->Internal )
680 std::string meshName;
681 vtkDataSet *ret(this->Internal->Tree.buildVTKInstance(this->Internal->IsStdOrMode,reqTS,meshName,this->Internal->TK,internalInfo));
682 if(this->Internal->GenerateVect)
684 vtkGenerateVectors::Operate(ret->GetPointData());
685 vtkGenerateVectors::Operate(ret->GetCellData());
686 vtkGenerateVectors::Operate(ret->GetFieldData());
687 // The operations above have potentially created new arrays -> This breaks the optimization of StaticMesh that expects the same field arrays over time.
688 // To enforce the cache recomputation declare modification of mesh.
689 //vtkGenerateVectors::ChangeMeshTimeToUpdateCache(ret);
694 void vtkMEDReader::PrintSelf(ostream& os, vtkIndent indent)
696 this->Superclass::PrintSelf(os, indent);