1 // Copyright (C) 2021 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 "vtkElectromagnetismRotation.h"
23 #include "MEDCouplingMemArray.hxx"
25 #include "ElectromagnetismRotationHelper.h"
26 #include "VTKMEDTraits.hxx"
28 #include "vtkAdjacentVertexIterator.h"
29 #include "vtkAOSDataArrayTemplate.h"
30 #include "vtkIntArray.h"
31 #include "vtkLongArray.h"
33 #include "vtkLongLongArray.h"
35 #include "vtkCellData.h"
36 #include "vtkPointData.h"
38 #include "vtkStreamingDemandDrivenPipeline.h"
39 #include "vtkUnstructuredGrid.h"
40 #include "vtkMultiBlockDataSet.h"
42 #include "vtkInformationStringKey.h"
43 #include "vtkAlgorithmOutput.h"
44 #include "vtkObjectFactory.h"
45 #include "vtkMutableDirectedGraph.h"
46 #include "vtkMultiBlockDataSet.h"
47 #include "vtkDataSet.h"
48 #include "vtkInformationVector.h"
49 #include "vtkInformation.h"
50 #include "vtkDataArraySelection.h"
51 #include "vtkTimeStamp.h"
52 #include "vtkInEdgeIterator.h"
53 #include "vtkInformationDataObjectKey.h"
54 #include "vtkExecutive.h"
55 #include "vtkVariantArray.h"
56 #include "vtkStringArray.h"
57 #include "vtkDoubleArray.h"
58 #include "vtkCharArray.h"
59 #include "vtkUnsignedCharArray.h"
60 #include "vtkDataSetAttributes.h"
61 #include "vtkDemandDrivenPipeline.h"
62 #include "vtkDataObjectTreeIterator.h"
63 #include "vtkThreshold.h"
64 #include "vtkMultiBlockDataGroupFilter.h"
65 #include "vtkCompositeDataToUnstructuredGridFilter.h"
66 #include "vtkInformationDataObjectMetaDataKey.h"
67 #include "vtkTransform.h"
68 #include "vtkFunctionParser.h"
75 const char ZE_SEP[]="@@][@@";
77 const char TS_STR[]="TS";
79 const char COM_SUP_STR[]="ComSup";
81 const char FAMILY_ID_CELL_NAME[]="FamilyIdCell";
83 const char NUM_ID_CELL_NAME[]="NumIdCell";
85 const char FAMILY_ID_NODE_NAME[]="FamilyIdNode";
87 const char NUM_ID_NODE_NAME[]="NumIdNode";
89 const char GLOBAL_NODE_ID_NAME[]="GlobalNodeIds";
91 vtkStandardNewMacro(vtkElectromagnetismRotation);
93 vtkInformationDataObjectMetaDataKey* GetMEDReaderMetaDataIfAny()
95 static const char ZE_KEY[] = "vtkMEDReader::META_DATA";
96 MEDCoupling::GlobalDict* gd(MEDCoupling::GlobalDict::GetInstance());
97 if (!gd->hasKey(ZE_KEY))
99 std::string ptSt(gd->value(ZE_KEY));
101 std::istringstream iss(ptSt);
103 return reinterpret_cast<vtkInformationDataObjectMetaDataKey*>(pt);
106 bool IsInformationOK(vtkInformation* info)
108 vtkInformationDataObjectMetaDataKey* key(GetMEDReaderMetaDataIfAny());
111 // Check the information contain meta data key
115 vtkMutableDirectedGraph* sil(vtkMutableDirectedGraph::SafeDownCast(info->Get(key)));
119 vtkAbstractArray* verticesNames(sil->GetVertexData()->GetAbstractArray("Names", idNames));
120 vtkStringArray* verticesNames2(vtkStringArray::SafeDownCast(verticesNames));
123 for (int i = 0; i < verticesNames2->GetNumberOfValues(); i++)
125 vtkStdString& st(verticesNames2->GetValue(i));
126 if (st == "MeshesFamsGrps")
132 class vtkElectromagnetismRotation::vtkElectromagnetismRotationInternal : public ElectromagnetismRotationInternal
138 vtkElectromagnetismRotation::vtkElectromagnetismRotation():SIL(NULL),Internal(new vtkElectromagnetismRotationInternal),InsideOut(0),Axis(2),RotationRotor(1e300)
142 vtkElectromagnetismRotation::~vtkElectromagnetismRotation()
144 delete this->Internal;
147 void vtkElectromagnetismRotation::SetInsideOut(int val)
149 if(this->InsideOut!=val)
156 void vtkElectromagnetismRotation::SetAxis(int axis)
165 void vtkElectromagnetismRotation::SetAngularStep(char *angStep)
167 this->AngularStep = angStep;
171 int vtkElectromagnetismRotation::RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
173 // vtkUnstructuredGridAlgorithm::RequestInformation(request,inputVector,outputVector);
176 // std::cerr << "########################################## vtkElectromagnetismRotation::RequestInformation ##########################################" << std::endl;
177 // request->Print(cout);
178 vtkInformation *outInfo(outputVector->GetInformationObject(0));
179 vtkInformation *inputInfo(inputVector[0]->GetInformationObject(0));
180 if(!ElectromagnetismRotationInternal::IndependantIsInformationOK(GetMEDReaderMetaDataIfAny(),inputInfo))
182 vtkErrorMacro("No SIL Data available ! The source of this filter must be MEDReader !");
186 this->SetSIL(vtkMutableDirectedGraph::SafeDownCast(inputInfo->Get(GetMEDReaderMetaDataIfAny())));
187 this->Internal->loadFrom(this->SIL);
188 //this->Internal->printMySelf(std::cerr);
190 catch(INTERP_KERNEL::Exception& e)
192 std::cerr << "Exception has been thrown in vtkElectromagnetismRotation::RequestInformation : " << e.what() << std::endl;
199 * Do not use vtkCxxSetObjectMacro macro because input mdg comes from an already managed in the pipeline just a ref on it.
201 void vtkElectromagnetismRotation::SetSIL(vtkMutableDirectedGraph *mdg)
208 template<class CellPointExtractor>
209 vtkDataSet *FilterFamilies(vtkThreshold *thres,
210 vtkDataSet *input, const std::set<int>& idsToKeep, bool insideOut, const char *arrNameOfFamilyField,
211 const char *associationForThreshold, bool& catchAll, bool& catchSmth)
213 const int VTK_DATA_ARRAY_DELETE=vtkAOSDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
214 const char ZE_SELECTION_ARR_NAME[]="@@ZeSelection@@";
215 vtkDataSet *output(input->NewInstance());
216 output->ShallowCopy(input);
217 thres->SetInputData(output);
218 vtkDataSetAttributes *dscIn(input->GetCellData()),*dscIn2(input->GetPointData());
219 vtkDataSetAttributes *dscOut(output->GetCellData()),*dscOut2(output->GetPointData());
221 double vMin(insideOut==0?1.:0.),vMax(insideOut==0?2.:1.);
222 thres->ThresholdBetween(vMin,vMax);
225 CellPointExtractor cpe2(input);
226 vtkDataArray *da(cpe2.Get()->GetScalars(arrNameOfFamilyField));
229 std::string daName(da->GetName());
230 typedef MEDFileVTKTraits<mcIdType>::VtkType vtkMCIdTypeArray;
231 vtkMCIdTypeArray *dai(vtkMCIdTypeArray::SafeDownCast(da));
232 if(daName!=arrNameOfFamilyField || !dai)
235 int nbOfTuples(dai->GetNumberOfTuples());
236 vtkCharArray *zeSelection(vtkCharArray::New());
237 zeSelection->SetName(ZE_SELECTION_ARR_NAME);
238 zeSelection->SetNumberOfComponents(1);
239 char *pt(new char[nbOfTuples]);
240 zeSelection->SetArray(pt,nbOfTuples,0,VTK_DATA_ARRAY_DELETE);
241 const mcIdType *inPtr(dai->GetPointer(0));
242 std::fill(pt,pt+nbOfTuples,0);
243 catchAll=true; catchSmth=false;
244 std::vector<bool> pt2(nbOfTuples,false);
245 for(std::set<int>::const_iterator it=idsToKeep.begin();it!=idsToKeep.end();it++)
247 bool catchFid(false);
248 for(int i=0;i<nbOfTuples;i++)
250 { pt2[i]=true; catchFid=true; }
256 for(int ii=0;ii<nbOfTuples;ii++)
259 CellPointExtractor cpe3(output);
260 int idx(cpe3.Get()->AddArray(zeSelection));
261 cpe3.Get()->SetActiveAttribute(idx,vtkDataSetAttributes::SCALARS);
262 cpe3.Get()->CopyScalarsOff();
263 zeSelection->Delete();
265 thres->SetInputArrayToProcess(idx,0,0,associationForThreshold,ZE_SELECTION_ARR_NAME);
267 vtkUnstructuredGrid *zeComputedOutput(thres->GetOutput());
268 CellPointExtractor cpe(zeComputedOutput);
269 cpe.Get()->RemoveArray(idx);
271 zeComputedOutput->Register(0);
272 return zeComputedOutput;
278 CellExtractor(vtkDataSet *ds):_ds(ds) { }
279 vtkDataSetAttributes *Get() { return _ds->GetCellData(); }
287 PointExtractor(vtkDataSet *ds):_ds(ds) { }
288 vtkDataSetAttributes *Get() { return _ds->GetPointData(); }
293 int vtkElectromagnetismRotation::RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
297 // std::cerr << "########################################## vtkElectromagnetismRotation::RequestData ##########################################" << std::endl;
298 // request->Print(cout);
299 vtkInformation* inputInfo=inputVector[0]->GetInformationObject(0);
300 vtkMultiBlockDataSet *inputMB(vtkMultiBlockDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
301 if(inputMB->GetNumberOfBlocks()!=1)
303 vtkErrorMacro(<< "vtkElectromagnetismRotation::RequestData : input has not the right number of parts ! Expected 1 !" ) ;
306 vtkDataSet *input(vtkDataSet::SafeDownCast(inputMB->GetBlock(0)));
307 vtkInformation *info(input->GetInformation());
308 vtkInformation *outInfo(outputVector->GetInformationObject(0));
309 vtkMultiBlockDataSet *output(vtkMultiBlockDataSet::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
310 output->SetNumberOfBlocks(2);
311 std::set<int> idsToKeep(this->Internal->getIdsToKeep());
312 this->Internal->clearSelection();
313 // first shrink the input
314 bool catchAll,catchSmth;
315 vtkNew<vtkThreshold> thres1,thres2;
316 vtkSmartPointer<vtkDataSet> rotor(FilterFamilies<CellExtractor>(thres1,input,idsToKeep,0,FAMILY_ID_CELL_NAME,"vtkDataObject::FIELD_ASSOCIATION_CELLS",catchAll,catchSmth));
317 vtkSmartPointer<vtkDataSet> stator(FilterFamilies<CellExtractor>(thres2,input,idsToKeep,1,FAMILY_ID_CELL_NAME,"vtkDataObject::FIELD_ASSOCIATION_CELLS",catchAll,catchSmth));
321 std::unique_ptr<double[]> timeSteps;
323 if(outInfo->Has(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP()))
324 reqTS=outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_TIME_STEP());
325 if(outInfo->Has(vtkStreamingDemandDrivenPipeline::TIME_STEPS()))
327 nbOfSteps = outInfo->Length(vtkStreamingDemandDrivenPipeline::TIME_STEPS());
328 timeSteps.reset(new double[ nbOfSteps ]);
329 outInfo->Get(vtkStreamingDemandDrivenPipeline::TIME_STEPS(),timeSteps.get());
330 //std::cerr << "nb : " << nbOfSteps << std::endl;
331 //std::for_each(timeSteps.get(),timeSteps.get()+nbOfSteps,[](double v) { std::cerr << v << std::endl; });
333 if(nbOfSteps<2 || !timeSteps.get())
335 vtkErrorMacro(<< "vtkElectromagnetismRotation::RequestData : A temporal dataset is expected ! Here < 2 time steps found !" ) ;
338 // Calcul de l angle effectif de rotation
339 double minTime(timeSteps[0]),maxTime(timeSteps[nbOfSteps-1]);
340 if(minTime == maxTime)
342 vtkErrorMacro(<< "vtkElectromagnetismRotation::RequestData : minTime == maxTime !" ) ;
345 double angularStepD = 0;
347 vtkNew<vtkFunctionParser> fp;
348 fp->SetFunction(this->AngularStep.c_str());
349 angularStepD = fp->GetScalarResult();
351 double effectiveTime(reqTS); effectiveTime = std::max(effectiveTime,minTime); effectiveTime = std::min(effectiveTime,maxTime);
352 this->RotationRotor = (angularStepD*((effectiveTime-minTime)/(maxTime-minTime)));
353 //std::cout << "*** " << effectiveTime << " " << minTime << " " << maxTime << " " << angleDegree << std::endl << std::flush;
354 //std::cout << "*** " << this->RotationRotor << std::endl << std::flush;
360 vtkNew<vtkTransform> transformR;
365 transformR->RotateX(this->RotationRotor);
370 transformR->RotateY(this->RotationRotor);
375 transformR->RotateZ(this->RotationRotor);
380 vtkErrorMacro(<< "vtkElectromagnetismRotation::RequestData : not recognized axis !" ) ;
384 vtkNew<vtkPoints> newCoords;
385 transformR->TransformPoints(vtkPointSet::SafeDownCast(rotor)->GetPoints(),newCoords);
386 vtkPointSet::SafeDownCast(rotor)->SetPoints(newCoords);
387 output->SetBlock(0,rotor);
388 output->SetBlock(1,stator);
395 catch(INTERP_KERNEL::Exception& e)
397 vtkErrorMacro(<< "Exception has been thrown in vtkElectromagnetismRotation::RequestData : " << e.what());
402 int vtkElectromagnetismRotation::GetSILUpdateStamp()
404 return (int)this->SILTime;
407 const char* vtkElectromagnetismRotation::GetGrpStart()
409 return ElectromagnetismRotationGrp::start();
412 const char* vtkElectromagnetismRotation::GetFamStart()
414 return ElectromagnetismRotationFam::start();
417 void vtkElectromagnetismRotation::PrintSelf(ostream& os, vtkIndent indent)
419 this->Superclass::PrintSelf(os, indent);
422 int vtkElectromagnetismRotation::GetNumberOfGroupsFlagsArrays()
424 int ret(this->Internal->getNumberOfEntries());
425 //std::cerr << "vtkElectromagnetismRotation::GetNumberOfFieldsTreeArrays() -> " << ret << std::endl;
429 const char *vtkElectromagnetismRotation::GetGroupsFlagsArrayName(int index)
431 const char *ret(this->Internal->getKeyOfEntry(index));
432 // std::cerr << "vtkElectromagnetismRotation::GetFieldsTreeArrayName(" << index << ") -> " << ret << std::endl;
436 int vtkElectromagnetismRotation::GetGroupsFlagsArrayStatus(const char *name)
438 int ret((int)this->Internal->getStatusOfEntryStr(name));
439 // std::cerr << "vtkElectromagnetismRotation::GetGroupsFlagsArrayStatus(" << name << ") -> " << ret << std::endl;
443 void vtkElectromagnetismRotation::SetGroupsFlagsStatus(const char *name, int status)
445 //std::cerr << "vtkElectromagnetismRotation::SetFieldsStatus(" << name << "," << status << ")" << std::endl;
446 this->Internal->setStatusOfEntryStr(name,(bool)status);
448 //this->Internal->printMySelf(std::cerr);
451 const char *vtkElectromagnetismRotation::GetMeshName()
453 return this->Internal->getMeshName();