1 // Copyright (C) 2010-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 "vtkExtractGroup.h"
22 #include "MEDFileFieldRepresentationTree.hxx"
23 #include "ExtractGroupHelper.h"
24 #include "vtkMEDReader.h"
25 #include "VTKMEDTraits.hxx"
27 #include "vtkAdjacentVertexIterator.h"
28 #include "vtkAOSDataArrayTemplate.h"
29 #include "vtkIntArray.h"
30 #include "vtkLongArray.h"
32 #include "vtkLongLongArray.h"
34 #include "vtkCellData.h"
35 #include "vtkPointData.h"
37 #include "vtkStreamingDemandDrivenPipeline.h"
38 #include "vtkUnstructuredGrid.h"
39 #include "vtkMultiBlockDataSet.h"
41 #include "vtkInformationStringKey.h"
42 #include "vtkAlgorithmOutput.h"
43 #include "vtkObjectFactory.h"
44 #include "vtkMutableDirectedGraph.h"
45 #include "vtkMultiBlockDataSet.h"
46 #include "vtkDataSet.h"
47 #include "vtkInformationVector.h"
48 #include "vtkInformation.h"
49 #include "vtkDataArraySelection.h"
50 #include "vtkTimeStamp.h"
51 #include "vtkInEdgeIterator.h"
52 #include "vtkInformationDataObjectKey.h"
53 #include "vtkExecutive.h"
54 #include "vtkVariantArray.h"
55 #include "vtkStringArray.h"
56 #include "vtkDoubleArray.h"
57 #include "vtkCharArray.h"
58 #include "vtkUnsignedCharArray.h"
59 #include "vtkDataSetAttributes.h"
60 #include "vtkDemandDrivenPipeline.h"
61 #include "vtkDataObjectTreeIterator.h"
62 #include "vtkThreshold.h"
63 #include "vtkMultiBlockDataGroupFilter.h"
64 #include "vtkCompositeDataToUnstructuredGridFilter.h"
65 #include "vtkInformationDataObjectMetaDataKey.h"
70 vtkStandardNewMacro(vtkExtractGroup)
72 class vtkExtractGroup::vtkExtractGroupInternal : public ExtractGroupInternal
78 vtkExtractGroup::vtkExtractGroup():SIL(NULL),Internal(new vtkExtractGroupInternal),InsideOut(0)
82 vtkExtractGroup::~vtkExtractGroup()
84 delete this->Internal;
87 void vtkExtractGroup::SetInsideOut(int val)
89 if(this->InsideOut!=val)
96 int vtkExtractGroup::RequestInformation(vtkInformation * /*request*/, vtkInformationVector **inputVector, vtkInformationVector */*outputVector*/)
98 // vtkUnstructuredGridAlgorithm::RequestInformation(request,inputVector,outputVector);
101 // std::cerr << "########################################## vtkExtractGroup::RequestInformation ##########################################" << std::endl;
102 // request->Print(cout);
103 //vtkInformation *outInfo(outputVector->GetInformationObject(0)); // todo: unused
104 vtkInformation *inputInfo(inputVector[0]->GetInformationObject(0));
105 if(!ExtractGroupInternal::IndependantIsInformationOK(vtkMEDReader::META_DATA(),inputInfo))
107 vtkErrorMacro("No SIL Data available ! The source of this filter must be MEDReader !");
111 this->SetSIL(vtkMutableDirectedGraph::SafeDownCast(inputInfo->Get(vtkMEDReader::META_DATA())));
112 this->Internal->loadFrom(this->SIL);
113 //this->Internal->printMySelf(std::cerr);
115 catch(INTERP_KERNEL::Exception& e)
117 std::cerr << "Exception has been thrown in vtkExtractGroup::RequestInformation : " << e.what() << std::endl;
124 * Do not use vtkCxxSetObjectMacro macro because input mdg comes from an already managed in the pipeline just a ref on it.
126 void vtkExtractGroup::SetSIL(vtkMutableDirectedGraph *mdg)
133 template<class CellPointExtractor>
134 vtkDataSet *FilterFamilies(vtkSmartPointer<vtkThreshold>& thres,
135 vtkDataSet *input, const std::set<int>& idsToKeep, bool insideOut, const char *arrNameOfFamilyField,
136 const char *associationForThreshold, bool& catchAll, bool& catchSmth)
138 const int VTK_DATA_ARRAY_DELETE=vtkAOSDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
139 const char ZE_SELECTION_ARR_NAME[]="@@ZeSelection@@";
140 vtkDataSet *output(input->NewInstance());
141 output->ShallowCopy(input);
142 thres->SetInputData(output);
143 //vtkDataSetAttributes *dscIn(input->GetCellData()),*dscIn2(input->GetPointData()); // todo: unused
144 //vtkDataSetAttributes *dscOut(output->GetCellData()),*dscOut2(output->GetPointData()); // todo: unused
146 double vMin(insideOut==0?1.:0.),vMax(insideOut==0?2.:1.);
147 thres->ThresholdBetween(vMin,vMax);
150 CellPointExtractor cpe2(input);
151 vtkDataArray *da(cpe2.Get()->GetScalars(arrNameOfFamilyField));
154 std::string daName(da->GetName());
155 typedef MEDFileVTKTraits<MEDCoupling::mcIdType>::VtkType vtkMCIdTypeArray;
156 vtkMCIdTypeArray *dai(vtkMCIdTypeArray::SafeDownCast(da));
157 if(daName!=arrNameOfFamilyField || !dai)
160 int nbOfTuples(dai->GetNumberOfTuples());
161 vtkCharArray *zeSelection(vtkCharArray::New());
162 zeSelection->SetName(ZE_SELECTION_ARR_NAME);
163 zeSelection->SetNumberOfComponents(1);
164 char *pt(new char[nbOfTuples]);
165 zeSelection->SetArray(pt,nbOfTuples,0,VTK_DATA_ARRAY_DELETE);
166 const MEDCoupling::mcIdType *inPtr(dai->GetPointer(0));
167 std::fill(pt,pt+nbOfTuples,0);
168 catchAll=true; catchSmth=false;
169 std::vector<bool> pt2(nbOfTuples,false);
170 for(std::set<int>::const_iterator it=idsToKeep.begin();it!=idsToKeep.end();it++)
172 bool catchFid(false);
173 for(int i=0;i<nbOfTuples;i++)
175 { pt2[i]=true; catchFid=true; }
181 for(int ii=0;ii<nbOfTuples;ii++)
184 CellPointExtractor cpe3(output);
185 int idx(cpe3.Get()->AddArray(zeSelection));
186 cpe3.Get()->SetActiveAttribute(idx,vtkDataSetAttributes::SCALARS);
187 cpe3.Get()->CopyScalarsOff();
188 zeSelection->Delete();
190 thres->SetInputArrayToProcess(idx,0,0,associationForThreshold,ZE_SELECTION_ARR_NAME);
192 vtkUnstructuredGrid *zeComputedOutput(thres->GetOutput());
193 CellPointExtractor cpe(zeComputedOutput);
194 cpe.Get()->RemoveArray(idx);
196 zeComputedOutput->Register(0);
197 return zeComputedOutput;
203 CellExtractor(vtkDataSet *ds):_ds(ds) { }
204 vtkDataSetAttributes *Get() { return _ds->GetCellData(); }
212 PointExtractor(vtkDataSet *ds):_ds(ds) { }
213 vtkDataSetAttributes *Get() { return _ds->GetPointData(); }
217 int vtkExtractGroup::RequestData(vtkInformation * /*request*/, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
221 // std::cerr << "########################################## vtkExtractGroup::RequestData ##########################################" << std::endl;
222 // request->Print(cout);
223 vtkInformation* inputInfo=inputVector[0]->GetInformationObject(0);
224 vtkMultiBlockDataSet *inputMB(vtkMultiBlockDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
225 if(inputMB->GetNumberOfBlocks()!=1)
227 std::ostringstream oss; oss << "vtkExtractGroup::RequestData : input has not the right number of parts ! Expected 1 !";
228 if(this->HasObserver("ErrorEvent") )
229 this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
231 vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
232 vtkObject::BreakOnError();
235 vtkDataSet *input(vtkDataSet::SafeDownCast(inputMB->GetBlock(0)));
236 //vtkInformation *info(input->GetInformation()); // todo: unused
237 vtkInformation *outInfo(outputVector->GetInformationObject(0));
238 vtkMultiBlockDataSet *output(vtkMultiBlockDataSet::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
239 std::set<int> idsToKeep(this->Internal->getIdsToKeep());
240 this->Internal->clearSelection();
241 // first shrink the input
242 bool catchAll,catchSmth;
243 vtkSmartPointer<vtkThreshold> thres1(vtkSmartPointer<vtkThreshold>::New()),thres2(vtkSmartPointer<vtkThreshold>::New());
244 vtkDataSet *tryOnCell(FilterFamilies<CellExtractor>(thres1,input,idsToKeep,this->InsideOut,
245 MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME,"vtkDataObject::FIELD_ASSOCIATION_CELLS",catchAll,catchSmth));
250 output->SetBlock(0,tryOnCell);
251 tryOnCell->Delete();//
258 vtkDataSet *tryOnNode(FilterFamilies<PointExtractor>(thres2,input,idsToKeep,this->InsideOut,
259 MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME,"vtkDataObject::FIELD_ASSOCIATION_POINTS",catchAll,catchSmth));
260 if(tryOnNode && catchSmth)
262 output->SetBlock(0,tryOnCell);
263 output->SetBlock(1,tryOnNode);
272 output->SetBlock(0,tryOnCell);
279 vtkDataSet *tryOnNode(FilterFamilies<PointExtractor>(thres1,input,idsToKeep,this->InsideOut,
280 MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME,"vtkDataObject::FIELD_ASSOCIATION_POINTS",catchAll,catchSmth));
284 output->SetBlock(0,tryOnNode);
290 output->SetBlock(0,tryOnNode);
299 vtkDataSet *tryOnNode(FilterFamilies<PointExtractor>(thres1,input,idsToKeep,this->InsideOut,
300 MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME,"vtkDataObject::FIELD_ASSOCIATION_POINTS",catchAll,catchSmth));
303 output->ShallowCopy(tryOnNode);
304 tryOnNode->Delete();//
309 std::ostringstream oss; oss << "vtkExtractGroup::RequestData : The integer array with name \""<< MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_CELL_NAME;
310 oss << "\" or \"" << MEDFileFieldRepresentationLeavesArrays::FAMILY_ID_NODE_NAME << "\" does not exist ! The extraction of group and/or family is not possible !";
311 if(this->HasObserver("ErrorEvent") )
312 this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
314 vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
315 vtkObject::BreakOnError();
320 catch(INTERP_KERNEL::Exception& e)
322 std::cerr << "Exception has been thrown in vtkExtractGroup::RequestData : " << e.what() << std::endl;
327 int vtkExtractGroup::GetSILUpdateStamp()
329 return (int)this->SILTime;
332 const char* vtkExtractGroup::GetGrpStart()
334 return ExtractGroupGrp::START;
337 const char* vtkExtractGroup::GetFamStart()
339 return ExtractGroupFam::START;
342 void vtkExtractGroup::PrintSelf(ostream& os, vtkIndent indent)
344 this->Superclass::PrintSelf(os, indent);
347 int vtkExtractGroup::GetNumberOfGroupsFlagsArrays()
349 int ret(this->Internal->getNumberOfEntries());
350 //std::cerr << "vtkExtractGroup::GetNumberOfFieldsTreeArrays() -> " << ret << std::endl;
354 const char *vtkExtractGroup::GetGroupsFlagsArrayName(int index)
356 const char *ret(this->Internal->getKeyOfEntry(index));
357 // std::cerr << "vtkExtractGroup::GetFieldsTreeArrayName(" << index << ") -> " << ret << std::endl;
361 int vtkExtractGroup::GetGroupsFlagsArrayStatus(const char *name)
363 int ret((int)this->Internal->getStatusOfEntryStr(name));
364 // std::cerr << "vtkExtractGroup::GetGroupsFlagsArrayStatus(" << name << ") -> " << ret << std::endl;
368 void vtkExtractGroup::SetGroupsFlagsStatus(const char *name, int status)
370 //std::cerr << "vtkExtractGroup::SetFieldsStatus(" << name << "," << status << ")" << std::endl;
371 this->Internal->setStatusOfEntryStr(name,(bool)status);
373 //this->Internal->printMySelf(std::cerr);
376 const char *vtkExtractGroup::GetMeshName()
378 return this->Internal->getMeshName();