1 // Copyright (C) 2016-2023 CEA, EDF
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 (EDF R&D)
21 #include "vtkMEDWriter.h"
22 #include "VTKToMEDMem.h"
24 #include <vtkAdjacentVertexIterator.h>
25 #include <vtkAlgorithmOutput.h>
26 #include <vtkCellArray.h>
27 #include <vtkCellData.h>
28 #include <vtkCharArray.h>
29 #include <vtkDataArraySelection.h>
30 #include <vtkDataObjectTreeIterator.h>
31 #include <vtkDataSet.h>
32 #include <vtkDataSetAttributes.h>
33 #include <vtkDemandDrivenPipeline.h>
34 #include <vtkDoubleArray.h>
35 #include <vtkExecutive.h>
36 #include <vtkFloatArray.h>
37 #include <vtkInEdgeIterator.h>
38 #include <vtkInformation.h>
39 #include <vtkInformationDataObjectKey.h>
40 #include <vtkInformationDataObjectMetaDataKey.h>
41 #include <vtkInformationStringKey.h>
42 #include <vtkInformationVector.h>
43 #include <vtkIntArray.h>
44 #include <vtkLongArray.h>
45 #include <vtkMultiBlockDataSet.h>
46 #include <vtkMutableDirectedGraph.h>
47 #include <vtkObjectFactory.h>
48 #include <vtkPointData.h>
49 #include <vtkPolyData.h>
50 #include <vtkRectilinearGrid.h>
51 #include <vtkStreamingDemandDrivenPipeline.h>
52 #include <vtkStringArray.h>
53 #include <vtkTimeStamp.h>
54 #include <vtkUnsignedCharArray.h>
55 #include <vtkUnstructuredGrid.h>
56 #include <vtkVariantArray.h>
57 #include <vtkWarpScalar.h>
58 #include <vtkAppendDataSets.h>
60 #include "MEDCouplingFieldDouble.hxx"
61 #include "MEDCouplingFieldFloat.hxx"
62 #include "MEDCouplingFieldInt.hxx"
63 #include "MEDCouplingMemArray.hxx"
64 #include "MEDCouplingRefCountObject.hxx"
65 #include "MEDFileData.hxx"
66 #include "MEDFileField.hxx"
67 #include "MEDFileMesh.hxx"
68 #include "MEDLoaderTraits.hxx"
74 vtkStandardNewMacro(vtkMEDWriter)
76 using MEDCoupling::MCAuto;
77 using MEDCoupling::MEDFileData;
78 using VTKToMEDMemWriter::Fam;
79 using VTKToMEDMemWriter::Grp;
81 static vtkInformationDataObjectMetaDataKey* GetMEDReaderMetaDataIfAny()
83 static const char ZE_KEY[] = "vtkFileSeriesGroupReader::META_DATA";
84 MEDCoupling::GlobalDict* gd(MEDCoupling::GlobalDict::GetInstance());
85 if (!gd->hasKey(ZE_KEY))
87 std::string ptSt(gd->value(ZE_KEY));
89 std::istringstream iss(ptSt);
91 return reinterpret_cast<vtkInformationDataObjectMetaDataKey*>(pt);
94 static bool IsInformationOK(vtkInformation* info)
96 vtkInformationDataObjectMetaDataKey* key(GetMEDReaderMetaDataIfAny());
99 // Check the information contain meta data key
103 vtkMutableDirectedGraph* sil(vtkMutableDirectedGraph::SafeDownCast(info->Get(key)));
107 vtkAbstractArray* verticesNames(sil->GetVertexData()->GetAbstractArray("Names", idNames));
108 vtkStringArray* verticesNames2(vtkStringArray::SafeDownCast(verticesNames));
111 for (int i = 0; i < verticesNames2->GetNumberOfValues(); i++)
113 vtkStdString& st(verticesNames2->GetValue(i));
114 if (st == "MeshesFamsGrps")
120 static void LoadFamGrpMapInfo(vtkMutableDirectedGraph* sil, std::string& meshName,
121 std::vector<Grp>& groups, std::vector<Fam>& fams)
124 throw MZCException("LoadFamGrpMapInfo : internal error !");
126 vtkAbstractArray* verticesNames(sil->GetVertexData()->GetAbstractArray("Names", idNames));
127 vtkStringArray* verticesNames2(vtkStringArray::SafeDownCast(verticesNames));
130 for (int i = 0; i < verticesNames2->GetNumberOfValues(); i++)
132 vtkStdString& st(verticesNames2->GetValue(i));
133 if (st == "MeshesFamsGrps")
140 throw INTERP_KERNEL::Exception(
141 "There is an internal error ! The tree on server side has not the expected look !");
142 vtkAdjacentVertexIterator* it0(vtkAdjacentVertexIterator::New());
143 sil->GetAdjacentVertices(id0, it0);
144 while (it0->HasNext())
146 vtkIdType id1(it0->Next());
147 std::string mName((const char*)verticesNames2->GetValue(id1));
149 vtkAdjacentVertexIterator* it1(vtkAdjacentVertexIterator::New());
150 sil->GetAdjacentVertices(id1, it1);
151 vtkIdType idZeGrps(it1->Next()); // zeGroups
152 vtkAdjacentVertexIterator* itGrps(vtkAdjacentVertexIterator::New());
153 sil->GetAdjacentVertices(idZeGrps, itGrps);
154 while (itGrps->HasNext())
156 vtkIdType idg(itGrps->Next());
157 Grp grp((const char*)verticesNames2->GetValue(idg));
158 vtkAdjacentVertexIterator* itGrps2(vtkAdjacentVertexIterator::New());
159 sil->GetAdjacentVertices(idg, itGrps2);
160 std::vector<std::string> famsOnGroup;
161 while (itGrps2->HasNext())
163 vtkIdType idgf(itGrps2->Next());
164 famsOnGroup.push_back(std::string((const char*)verticesNames2->GetValue(idgf)));
166 grp.setFamilies(famsOnGroup);
168 groups.push_back(grp);
171 vtkIdType idZeFams(it1->Next()); // zeFams
173 vtkAdjacentVertexIterator* itFams(vtkAdjacentVertexIterator::New());
174 sil->GetAdjacentVertices(idZeFams, itFams);
175 while (itFams->HasNext())
177 vtkIdType idf(itFams->Next());
178 Fam fam((const char*)verticesNames2->GetValue(idf));
186 vtkMEDWriter::vtkMEDWriter()
187 : UnPolygonize(0),WriteAllTimeSteps(0)
188 , NumberOfTimeSteps(0)
189 , CurrentTimeIndex(0)
194 vtkMEDWriter::~vtkMEDWriter() {}
196 int vtkMEDWriter::RequestInformation(
197 vtkInformation* /*request*/, vtkInformationVector** inputVector, vtkInformationVector* /*outputVector*/)
199 // std::cerr << "########################################## vtkMEDWriter::RequestInformation
200 // ########################################## " << (const char *) this->FileName << std::endl;
201 vtkInformation* inInfo(inputVector[0]->GetInformationObject(0));
202 if (inInfo->Has(vtkStreamingDemandDrivenPipeline::TIME_STEPS()))
203 this->NumberOfTimeSteps = inInfo->Length(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
205 this->NumberOfTimeSteps = 0;
209 int vtkMEDWriter::RequestUpdateExtent(vtkInformation* vtkNotUsed(request),
210 vtkInformationVector** inputVector, vtkInformationVector* vtkNotUsed(outputVector))
213 inputVector[0]->GetInformationObject(0)->Get(vtkStreamingDemandDrivenPipeline::TIME_STEPS()));
214 if (inTimes && this->WriteAllTimeSteps)
216 double timeReq(inTimes[this->CurrentTimeIndex]);
217 inputVector[0]->GetInformationObject(0)->Set(
218 vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP(), timeReq);
223 void ExceptionDisplayer(vtkMEDWriter* self, const std::string& fileName, std::exception& e)
225 std::ostringstream oss;
226 oss << "Exception has been thrown in vtkMEDWriter::RequestData : During writing of \"" << fileName
227 << "\", the following exception has been thrown : " << e.what() << std::endl;
228 if (self->HasObserver("ErrorEvent"))
229 self->InvokeEvent("ErrorEvent", const_cast<char*>(oss.str().c_str()));
231 vtkOutputWindowDisplayErrorText(const_cast<char*>(oss.str().c_str()));
232 vtkObject::BreakOnError();
235 int vtkMEDWriter::RequestData(
236 vtkInformation* request, vtkInformationVector** inputVector, vtkInformationVector* outputVector)
238 // std::cerr << "########################################## vtkMEDWriter::RequestData
239 // ########################################## " << (const char *) this->FileName << std::endl;
242 bool writeAll((this->WriteAllTimeSteps != 0 && this->NumberOfTimeSteps > 0));
245 if (this->CurrentTimeIndex == 0)
246 request->Set(vtkStreamingDemandDrivenPipeline::CONTINUE_EXECUTING(), 1);
250 request->Remove(vtkStreamingDemandDrivenPipeline::CONTINUE_EXECUTING());
251 this->CurrentTimeIndex = 0;
254 vtkInformation* inputInfo(inputVector[0]->GetInformationObject(0));
255 std::string meshName;
256 std::vector<Grp> groups;
257 std::vector<Fam> fams;
258 if (IsInformationOK(inputInfo))
260 vtkMutableDirectedGraph* famGrpGraph(
261 vtkMutableDirectedGraph::SafeDownCast(inputInfo->Get(GetMEDReaderMetaDataIfAny())));
262 LoadFamGrpMapInfo(famGrpGraph, meshName, groups, fams);
264 vtkInformation* outInfo(outputVector->GetInformationObject(0));
265 vtkDataObject* inputBase(vtkDataObject::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
268 "Not recognized data object in input of the MEDWriter ! Maybe not implemented yet !");
271 vtkInformation* inInfo(inputVector[0]->GetInformationObject(0));
272 vtkDataObject* input(vtkDataObject::GetData(inInfo));
273 timeStep = input->GetInformation()->Get(vtkDataObject::DATA_TIME_STEP());
276 vtkDataObject* input = inputBase;
277 vtkSmartPointer<vtkDataObject> inputUnpoly;
278 if(this->UnPolygonize)
280 vtkPointSet *inputPS = vtkPointSet::SafeDownCast(input);
282 throw MZCException("UnPolygonize is activated whereas it is not a PointSet type !");
283 vtkSmartPointer<vtkAppendDataSets> ads = vtkSmartPointer<vtkAppendDataSets>::New();
284 ads->SetInputData(inputPS);
286 vtkDataObject *inputPostPro = ads->GetOutputDataObject(0);
287 inputPostPro->Register(nullptr);
288 inputUnpoly.TakeReference(inputPostPro);
289 input = inputPostPro;
292 MCAuto<MEDFileData> mfd(MEDFileData::New());
293 VTKToMEDMemWriter::WriteMEDFileFromVTKGDS(mfd, input, timeStep, this->CurrentTimeIndex);
294 VTKToMEDMemWriter::PutFamGrpInfoIfAny(mfd, meshName, groups, fams);
295 if (this->WriteAllTimeSteps == 0 ||
296 (this->WriteAllTimeSteps != 0 && this->CurrentTimeIndex == 0))
297 mfd->write(this->FileName, 2);
300 mfd->getFields()->write(this->FileName, 0);
302 outInfo->Set(vtkDataObject::DATA_OBJECT(), input);
306 this->CurrentTimeIndex++;
307 if (this->CurrentTimeIndex >= this->NumberOfTimeSteps)
309 request->Remove(vtkStreamingDemandDrivenPipeline::CONTINUE_EXECUTING());
310 this->CurrentTimeIndex = 0;
314 catch (INTERP_KERNEL::Exception& e)
316 ExceptionDisplayer(this, (const char*)this->FileName, e);
319 catch (MZCException& e)
321 ExceptionDisplayer(this, (const char*)this->FileName, e);
327 void vtkMEDWriter::PrintSelf(ostream& os, vtkIndent indent)
329 this->Superclass::PrintSelf(os, indent);
332 int vtkMEDWriter::Write()