1 // Copyright (C) 2016-2022 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
21 #include "MEDFileField.hxx"
22 #include "MEDCouplingUMesh.hxx"
23 #include "MEDCouplingCMesh.hxx"
24 #include "MEDCouplingCurveLinearMesh.hxx"
25 #include "MEDFileMesh.hxx"
26 #include "MEDFileFieldRepresentationTree.hxx"
30 #include <vtkCPDataDescription.h>
31 #include <vtkCPInputDataDescription.h>
32 #include <vtkCPPythonScriptPipeline.h>
33 #include <vtkCPProcessor.h>
34 #include "vtkDataSet.h"
35 #include "vtkCellData.h"
37 #include <vtkSmartPointer.h>
38 #include <vtkXMLUnstructuredGridWriter.h>
39 #include <vtkUnstructuredGrid.h>
41 using namespace MEDCoupling;
44 void Visualization::CatalystInitialize(const std::string& script)
49 Processor = vtkCPProcessor::New();
50 Processor->Initialize();
54 Processor->RemoveAllPipelines();
56 // Add in the Python script
57 vtkCPPythonScriptPipeline* pipeline = vtkCPPythonScriptPipeline::New();
58 char *c = const_cast<char*>(script.c_str());
59 pipeline->Initialize(c);
60 Processor->AddPipeline(pipeline);
65 void Visualization::CatalystFinalize()
76 void Visualization::CatalystCoProcess(vtkDataSet *VTKGrid, double time,
77 unsigned int timeStep)
79 vtkCPDataDescription* dataDescription = vtkCPDataDescription::New();
80 // specify the simulation time and time step for Catalyst
81 dataDescription->AddInput("input");
82 dataDescription->SetTimeData(time, timeStep);
84 if(Processor->RequestDataDescription(dataDescription) != 0)
86 // Make a map from input to our VTK grid so that
87 // Catalyst gets the proper input dataset for the pipeline.
88 dataDescription->GetInputDescriptionByName("input")->SetGrid(VTKGrid);
89 // Call Catalyst to execute the desired pipelines.
90 Processor->CoProcess(dataDescription);
92 dataDescription->Delete();
95 void Visualization::ConvertToVTK(MEDCoupling::MEDCouplingFieldDouble* field, vtkDataSet *&VTKGrid)
96 //vtkDataSet * Visualization::ConvertToVTK(MEDCoupling::MEDCouplingFieldDouble* field)
98 MEDCoupling::MEDCouplingFieldDouble *f = field;
99 MEDCoupling::MCAuto<MEDFileField1TS> ff(MEDFileField1TS::New());
100 ff->setFieldNoProfileSBT(f);
101 MCAuto<MEDFileFieldMultiTS> fmts(MEDFileFieldMultiTS::New());
102 fmts->pushBackTimeStep(ff);
103 MCAuto<MEDFileFields> fs(MEDFileFields::New());
106 MEDCouplingMesh *m(f->getMesh());
107 MCAuto<MEDFileMesh> mm;
108 if(dynamic_cast<MEDCouplingUMesh *>(m))
110 MCAuto<MEDFileUMesh> mmu(MEDFileUMesh::New());
111 mmu->setMeshAtLevel(0,dynamic_cast<MEDCouplingUMesh *>(m));
112 mmu->forceComputationOfParts();
113 mm=DynamicCast<MEDFileUMesh,MEDFileMesh>(mmu);
115 else if(dynamic_cast<MEDCouplingCMesh *>(m))
117 MCAuto<MEDFileCMesh> mmc(MEDFileCMesh::New());
118 mmc->setMesh(dynamic_cast<MEDCouplingCMesh *>(m));
119 mm=DynamicCast<MEDFileCMesh,MEDFileMesh>(mmc);
121 else if(dynamic_cast<MEDCouplingCurveLinearMesh *>(m))
123 // MCAuto<MEDFileCLMesh> mmc(MEDFileCLMesh::New());
124 MCAuto<MEDFileCurveLinearMesh> mmc(MEDFileCurveLinearMesh::New());
125 mmc->setMesh(dynamic_cast<MEDCouplingCurveLinearMesh *>(m));
126 //mm=DynamicCast<MEDCouplingCurveLinearMesh,MEDFileMesh>(mmc);
133 MCAuto<MEDFileMeshes> ms(MEDFileMeshes::New());
135 ms->getMeshesNames();
138 MPI_Comm_rank(MPI_COMM_WORLD,&proc_id);
140 MEDFileFieldRepresentationTree Tree;
141 Tree.loadInMemory(fs,ms);
143 Tree.changeStatusOfAndUpdateToHaveCoherentVTKDataSet(0,true);
144 Tree.changeStatusOfAndUpdateToHaveCoherentVTKDataSet(1,false);
145 Tree.activateTheFirst();//"TS0/Fluid domain/ComSup0/TempC@@][@@P0"
146 std::string meshName;
149 double tmp3(f->getTime(tmp1,tmp2));
150 vtkDataSet *ret(Tree.buildVTKInstance(false,tmp3,meshName,TK));
154 Visualization::Visualization()
160 // ajouter en parametre le fichier python qui contient le pipeline
161 // enlever le const s'il gene
162 void Visualization::run(MEDCoupling::MEDCouplingFieldDouble* field, const std::string& pathPipeline)
165 MPI_Comm_rank(MPI_COMM_WORLD,&proc_id);
167 vtkDataSet *VTKGrid = 0;
168 ConvertToVTK(field, VTKGrid);
170 const char *fileName = pathPipeline.c_str();
171 CatalystInitialize(fileName);
172 CatalystCoProcess(VTKGrid, 0.0, 0);