src/Plugins/MEDReader/plugin/MEDReaderIO/vtkExtractCellType.cxx

   1 // Copyright (C) 2010-2021  CEA/DEN, EDF R&D
   2 //
   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.
   7 //
   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.
  12 //
  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
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay
  20
  21 #include "vtkExtractCellType.h"
  22 #include "MEDFileFieldRepresentationTree.hxx"
  23 #include "MEDFileFieldOverView.hxx"
  24
  25 #include "vtkAdjacentVertexIterator.h"
  26 #include "vtkAOSDataArrayTemplate.h"
  27 #include "vtkIntArray.h"
  28 #include "vtkCellData.h"
  29 #include "vtkPointData.h"
  30
  31 #include "vtkStreamingDemandDrivenPipeline.h"
  32 #include "vtkUnstructuredGrid.h"
  33 #include  "vtkMultiBlockDataSet.h"
  34
  35 #include "vtkInformationStringKey.h"
  36 #include "vtkAlgorithmOutput.h"
  37 #include "vtkObjectFactory.h"
  38 #include "vtkMutableDirectedGraph.h"
  39 #include "vtkMultiBlockDataSet.h"
  40 #include "vtkDataSet.h"
  41 #include "vtkInformationVector.h"
  42 #include "vtkInformation.h"
  43 #include "vtkDataArraySelection.h"
  44 #include "vtkTimeStamp.h"
  45 #include "vtkInEdgeIterator.h"
  46 #include "vtkInformationDataObjectKey.h"
  47 #include "vtkExecutive.h"
  48 #include "vtkVariantArray.h"
  49 #include "vtkStringArray.h"
  50 #include "vtkDoubleArray.h"
  51 #include "vtkCharArray.h"
  52 #include "vtkUnsignedCharArray.h"
  53 #include "vtkDataSetAttributes.h"
  54 #include "vtkDemandDrivenPipeline.h"
  55 #include "vtkDataObjectTreeIterator.h"
  56 #include "vtkThreshold.h"
  57
  58 #include <map>
  59 #include <deque>
  60
  61 vtkStandardNewMacro(vtkExtractCellType)
  62
  63 vtkCxxSetObjectMacro(vtkExtractCellType, SIL, vtkMutableDirectedGraph)
  64
  65 ///////////////////
  66
  67 /// \cond PRIVATE
  68 class ExtractCellTypeStatus
  69 {
  70 public:
  71   ExtractCellTypeStatus():_status(false),_vtkt(-1),_mct(INTERP_KERNEL::NORM_ERROR) { }
  72   ExtractCellTypeStatus(int vtkt, INTERP_KERNEL::NormalizedCellType mct);
  73   bool isSame(int vtkt, INTERP_KERNEL::NormalizedCellType mct) const { return _vtkt==vtkt && _mct==mct; }
  74   bool getStatus() const { return _status; }
  75   void setStatus(bool status) const { _status=status; }
  76   void cpyStatusFrom(const ExtractCellTypeStatus& other) { _status=other._status; }
  77   std::string getKey() const { return _type_str; }
  78   const char *getKeyOfEntry() const { return _type_str.c_str(); }
  79   int getVTKCellType() const { return _vtkt; }
  80   void printMySelf(std::ostream& os) const;
  81   bool isSameAs(const ExtractCellTypeStatus& other) const;
  82   void feedSIL(vtkMutableDirectedGraph *sil, vtkIdType root, vtkVariantArray *childEdge, std::vector<std::string>& names) const;
  83 protected:
  84   mutable bool _status;
  85   int _vtkt;
  86   INTERP_KERNEL::NormalizedCellType _mct;
  87   std::string _type_str;
  88 };
  89
  90 class vtkExtractCellType::vtkExtractCellTypeInternal
  91 {
  92 public:
  93   vtkExtractCellTypeInternal():_ref_mtime(0) { }
  94   int getNumberOfEntries() const;
  95   const char *getKeyOfEntry(int i) const;
  96   bool getStatusOfEntryStr(const char *entry) const;
  97   void setStatusOfEntryStr(const char *entry, bool status) const;
  98   void feedSIL(vtkMutableDirectedGraph *sil) const;
  99   std::vector<int> getIdsToKeep() const;
 100   void printMySelf(std::ostream& os) const;
 101   bool setRefTime(vtkObject *input) const;
 102   // non const methods
 103   void loadFrom(const std::map<int,INTERP_KERNEL::NormalizedCellType>& m);
 104 private:
 105   const ExtractCellTypeStatus& getEntry(const char *entry) const;
 106   bool checkSame(const std::map<int,INTERP_KERNEL::NormalizedCellType>& m) const;
 107 private:
 108   std::vector<ExtractCellTypeStatus> _types;
 109   mutable unsigned long _ref_mtime;
 110 };
 111
 112 bool vtkExtractCellType::vtkExtractCellTypeInternal::setRefTime(vtkObject *input) const
 113 {
 114   unsigned long mtime(input->GetMTime());
 115   if(mtime>_ref_mtime)
 116     {
 117       _ref_mtime=mtime;
 118       return true;
 119     }
 120   else
 121     return false;
 122 }
 123
 124 std::vector<int> vtkExtractCellType::vtkExtractCellTypeInternal::getIdsToKeep() const
 125 {
 126   std::vector<int> ret;
 127   for(std::vector<ExtractCellTypeStatus>::const_iterator it=_types.begin();it!=_types.end();it++)
 128     {
 129       if((*it).getStatus())
 130         ret.push_back((*it).getVTKCellType());
 131     }
 132   return ret;
 133 }
 134
 135 void vtkExtractCellType::vtkExtractCellTypeInternal::feedSIL(vtkMutableDirectedGraph *sil) const
 136 {
 137   vtkSmartPointer<vtkVariantArray> childEdge(vtkSmartPointer<vtkVariantArray>::New());
 138   childEdge->InsertNextValue(0);
 139   vtkSmartPointer<vtkVariantArray> crossEdge(vtkSmartPointer<vtkVariantArray>::New());
 140   crossEdge->InsertNextValue(1);
 141   // CrossEdge is an edge linking hierarchies.
 142   vtkUnsignedCharArray* crossEdgesArray=vtkUnsignedCharArray::New();
 143   crossEdgesArray->SetName("CrossEdges");
 144   sil->GetEdgeData()->AddArray(crossEdgesArray);
 145   crossEdgesArray->Delete();
 146   std::vector<std::string> names;
 147   // Add global fields root
 148   vtkIdType root(sil->AddVertex());
 149   names.push_back("CellTypesTree");
 150   //
 151   for(std::vector<ExtractCellTypeStatus>::const_iterator it=_types.begin();it!=_types.end();it++)
 152     {
 153       (*it).feedSIL(sil,root,childEdge,names);
 154     }
 155   // This array is used to assign names to nodes.
 156   vtkStringArray *namesArray(vtkStringArray::New());
 157   namesArray->SetName("Names");
 158   namesArray->SetNumberOfTuples(sil->GetNumberOfVertices());
 159   sil->GetVertexData()->AddArray(namesArray);
 160   namesArray->Delete();
 161   std::vector<std::string>::const_iterator iter;
 162   vtkIdType cc;
 163   for(cc=0, iter=names.begin(); iter!=names.end(); ++iter, ++cc)
 164     namesArray->SetValue(cc,(*iter).c_str());
 165 }
 166
 167 void vtkExtractCellType::vtkExtractCellTypeInternal::loadFrom(const std::map<int,INTERP_KERNEL::NormalizedCellType>& m)
 168 {
 169   if(checkSame(m))
 170     return;
 171   //
 172   std::size_t sz(m.size()),ii(0);
 173   _types.resize(sz);
 174   for(std::map<int,INTERP_KERNEL::NormalizedCellType>::const_iterator it=m.begin();it!=m.end();it++,ii++)
 175     {
 176       ExtractCellTypeStatus elt((*it).first,(*it).second);
 177       _types[ii]=elt;
 178     }
 179 }
 180
 181 int vtkExtractCellType::vtkExtractCellTypeInternal::getNumberOfEntries() const
 182 {
 183   return (int) _types.size();
 184 }
 185
 186 const char *vtkExtractCellType::vtkExtractCellTypeInternal::getKeyOfEntry(int i) const
 187 {
 188   return _types[i].getKeyOfEntry();
 189 }
 190
 191 bool vtkExtractCellType::vtkExtractCellTypeInternal::checkSame(const std::map<int,INTERP_KERNEL::NormalizedCellType>& m) const
 192 {
 193   std::size_t sz(m.size());
 194   if(sz!=_types.size())
 195     return false;
 196   bool ret(true);
 197   std::map<int,INTERP_KERNEL::NormalizedCellType>::const_iterator it(m.begin());
 198   for(std::size_t i=0;i<sz && ret;i++,it++)
 199     ret=_types[i].isSame((*it).first,(*it).second);
 200   return ret;
 201 }
 202
 203 const ExtractCellTypeStatus& vtkExtractCellType::vtkExtractCellTypeInternal::getEntry(const char *entry) const
 204 {
 205   std::string entryCpp(entry);
 206   for(std::vector<ExtractCellTypeStatus>::const_iterator it0=_types.begin();it0!=_types.end();it0++)
 207     if(entryCpp==(*it0).getKey())
 208       return *it0;
 209   std::ostringstream oss; oss << "vtkExtractCellTypeInternal::getEntry : no such entry \"" << entry << "\"!";
 210   throw INTERP_KERNEL::Exception(oss.str().c_str());
 211 }
 212
 213 bool vtkExtractCellType::vtkExtractCellTypeInternal::getStatusOfEntryStr(const char *entry) const
 214 {
 215   try
 216     {
 217       const ExtractCellTypeStatus& elt(getEntry(entry));
 218       return elt.getStatus();
 219     }
 220   catch (INTERP_KERNEL::Exception& /*e*/)
 221     {
 222       //std::cerr << vtkDebugMacro"Exception has been thrown in vtkExtractCellType::vtkExtractCellTypeInternal::getStatusOfEntryStr : " << e.what() << std::endl;
 223       return false;
 224     }
 225 }
 226
 227 void vtkExtractCellType::vtkExtractCellTypeInternal::setStatusOfEntryStr(const char *entry, bool status) const
 228 {
 229   try
 230     {
 231       const ExtractCellTypeStatus& elt(getEntry(entry));
 232       elt.setStatus(status);
 233     }
 234   catch (INTERP_KERNEL::Exception& /*e*/)
 235     {
 236       //std::cerr << "Exception has been thrown in vtkExtractCellType::vtkExtractCellTypeInternal::setStatusOfEntryStr : " << e.what() << std::endl;
 237     }
 238 }
 239
 240 void vtkExtractCellType::vtkExtractCellTypeInternal::printMySelf(std::ostream& os) const
 241 {
 242   for(std::vector<ExtractCellTypeStatus>::const_iterator it0=_types.begin();it0!=_types.end();it0++)
 243     (*it0).printMySelf(os);
 244 }
 245
 246 ExtractCellTypeStatus::ExtractCellTypeStatus(int vtkt, INTERP_KERNEL::NormalizedCellType mct):_status(false),_vtkt(vtkt),_mct(mct)
 247 {
 248   std::string name(INTERP_KERNEL::CellModel::GetCellModel(mct).getRepr());
 249   _type_str=name.substr(5);//skip "NORM_"
 250 }
 251
 252 void ExtractCellTypeStatus::printMySelf(std::ostream& os) const
 253 {
 254   os << "      -" << _type_str << "(";
 255   if(_status)
 256     os << "X";
 257   else
 258     os << " ";
 259   os << ")" << std::endl;
 260 }
 261
 262 bool ExtractCellTypeStatus::isSameAs(const ExtractCellTypeStatus& other) const
 263 {
 264   return _vtkt==other._vtkt && _mct==other._mct;
 265 }
 266
 267 void ExtractCellTypeStatus::feedSIL(vtkMutableDirectedGraph *sil, vtkIdType root, vtkVariantArray *childEdge, std::vector<std::string>& names) const
 268 {
 269   vtkIdType InfoGeoType(sil->AddChild(root,childEdge));
 270   names.push_back(_type_str);
 271   /*vtkIdType InfoVTKID(*/sil->AddChild(InfoGeoType,childEdge)/*)*/; // todo: unused
 272   std::ostringstream oss; oss << _vtkt;
 273   names.push_back(oss.str());
 274 }
 275 /// \endcond PRIVATE
 276
 277 ////////////////////
 278
 279 vtkExtractCellType::vtkExtractCellType():SIL(NULL),Internal(new vtkExtractCellTypeInternal),InsideOut(0)
 280 {
 281 }
 282
 283 vtkExtractCellType::~vtkExtractCellType()
 284 {
 285   if(this->SIL)
 286     this->SIL->Delete();
 287   delete this->Internal;
 288 }
 289
 290 void vtkExtractCellType::SetInsideOut(int val)
 291 {
 292   if(this->InsideOut!=val)
 293     {
 294       this->InsideOut=val;
 295       this->Modified();
 296     }
 297 }
 298
 299 int vtkExtractCellType::RequestInformation(vtkInformation * /*request*/, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 300 {
 301   try
 302     {
 303       //std::cerr << "########################################## vtkExtractCellType::RequestInformation ##########################################" << std::endl;
 304       vtkInformation *outInfo(outputVector->GetInformationObject(0));
 305       vtkInformation *inputInfo(inputVector[0]->GetInformationObject(0));
 306       vtkDataSet *input(0);
 307       {
 308         vtkDataObject *inp(inputInfo->Get(vtkDataObject::DATA_OBJECT()));
 309         if(vtkDataSet::SafeDownCast(inp))
 310           input=vtkDataSet::SafeDownCast(inp);
 311         else
 312           {
 313             vtkMultiBlockDataSet *inputTmp(vtkMultiBlockDataSet::SafeDownCast(inp));
 314             if(inputTmp)
 315               {
 316                 if(inputTmp->GetNumberOfBlocks()!=1)
 317                   {
 318                     vtkDebugMacro("vtkExtractCellType::RequestInformation : input vtkMultiBlockDataSet must contain exactly 1 block !");
 319                     return 0;
 320                   }
 321                 vtkDataSet *blk0(vtkDataSet::SafeDownCast(inputTmp->GetBlock(0)));
 322                 if(!blk0)
 323                   {
 324                     vtkDebugMacro("vtkExtractCellType::RequestInformation : the single block in input vtkMultiBlockDataSet must be a vtkDataSet instance !");
 325                     return 0;
 326                   }
 327                 input=blk0;
 328               }
 329             else
 330               {
 331                 vtkDebugMacro("vtkExtractCellType::RequestInformation : supported input are vtkDataSet or vtkMultiBlockDataSet !");
 332                 return 0;
 333               }
 334           }
 335       }
 336       if(this->Internal->setRefTime(input))
 337         {
 338           vtkIdType nbOfCells(input->GetNumberOfCells());
 339           std::map<int,INTERP_KERNEL::NormalizedCellType> m;
 340           for(vtkIdType cellId=0;cellId<nbOfCells;cellId++)
 341             {
 342               int vtkCt(input->GetCellType(cellId));
 343               const std::map<int,INTERP_KERNEL::NormalizedCellType>::const_iterator it(m.find(vtkCt));
 344               if(it==m.end())
 345                 {
 346                   const unsigned char *pos(std::find(MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH,vtkCt));
 347                   if(pos==MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE+MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE_LGTH)
 348                     {
 349                       vtkDebugMacro("vtkExtractCellType::RequestInformation : cell #" << cellId << " has unrecognized type !");
 350                       return 0;
 351                     }
 352                   m[vtkCt]=(INTERP_KERNEL::NormalizedCellType)std::distance(MEDCoupling::MEDMeshMultiLev::PARAMEDMEM_2_VTKTYPE,pos);
 353             }
 354             }
 355           this->Internal->loadFrom(m);
 356           if(this->SIL)
 357             this->SIL->Delete();
 358           this->SIL=vtkMutableDirectedGraph::New();
 359           this->Internal->feedSIL(this->SIL);
 360           //
 361           outInfo->Set(vtkDataObject::SIL(),this->SIL);
 362         }
 363     }
 364   catch(INTERP_KERNEL::Exception& e)
 365     {
 366       std::cerr << "Exception has been thrown in vtkExtractCellType::RequestInformation : " << e.what() << std::endl;
 367       return 0;
 368     }
 369   return 1;
 370 }
 371
 372 vtkDataSet *FilterFamilies(vtkDataSet *input, const std::vector<int>& idsToKeep, bool insideOut)
 373 {
 374   const int VTK_DATA_ARRAY_DELETE=vtkAOSDataArrayTemplate<double>::VTK_DATA_ARRAY_DELETE;
 375   const char ZE_SELECTION_ARR_NAME[]="@@ZeSelection@@";
 376   vtkDataSet *output(input->NewInstance());
 377   output->ShallowCopy(input);
 378   vtkSmartPointer<vtkThreshold> thres(vtkSmartPointer<vtkThreshold>::New());
 379   thres->SetInputData(output);
 380   //vtkDataSetAttributes *dscIn(input->GetCellData()),*dscIn2(input->GetPointData()); // todo: unused
 381   //vtkDataSetAttributes *dscOut(output->GetCellData()),*dscOut2(output->GetPointData()); // todo: unused
 382   //
 383   double vMin(insideOut==0?1.:0.),vMax(insideOut==0?2.:1.);
 384   thres->ThresholdBetween(vMin,vMax);
 385   // OK for the output
 386   vtkIdType nbOfCells(input->GetNumberOfCells());
 387   vtkCharArray *zeSelection(vtkCharArray::New());
 388   zeSelection->SetName(ZE_SELECTION_ARR_NAME);
 389   zeSelection->SetNumberOfComponents(1);
 390   char *pt(new char[nbOfCells]);
 391   zeSelection->SetArray(pt,nbOfCells,0,VTK_DATA_ARRAY_DELETE);
 392   std::fill(pt,pt+nbOfCells,0);
 393   std::vector<bool> pt2(nbOfCells,false);
 394   for(std::vector<int>::const_iterator it=idsToKeep.begin();it!=idsToKeep.end();it++)
 395     {
 396       for(vtkIdType ii=0;ii<nbOfCells;ii++)
 397         {
 398           if(input->GetCellType(ii)==*it)
 399             pt2[ii]=true;
 400         }
 401     }
 402   for(int ii=0;ii<nbOfCells;ii++)
 403     if(pt2[ii])
 404       pt[ii]=2;
 405   int idx(output->GetCellData()->AddArray(zeSelection));
 406   output->GetCellData()->SetActiveAttribute(idx,vtkDataSetAttributes::SCALARS);
 407   output->GetCellData()->CopyScalarsOff();
 408   zeSelection->Delete();
 409   //
 410   thres->SetInputArrayToProcess(idx,0,0,"vtkDataObject::FIELD_ASSOCIATION_CELLS",ZE_SELECTION_ARR_NAME);
 411   thres->Update();
 412   vtkUnstructuredGrid *zeComputedOutput(thres->GetOutput());
 413   zeComputedOutput->GetCellData()->RemoveArray(idx);
 414   output->Delete();
 415   zeComputedOutput->Register(0);
 416   return zeComputedOutput;
 417 }
 418
 419 int vtkExtractCellType::RequestData(vtkInformation * /*request*/, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 420 {
 421   try
 422     {
 423       //std::cerr << "########################################## vtkExtractCellType::RequestData        ##########################################" << std::endl;
 424       vtkInformation* inputInfo=inputVector[0]->GetInformationObject(0);
 425       vtkDataSet *input(vtkDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
 426       //vtkInformation *info(input->GetInformation()); // todo: unused
 427       vtkInformation *outInfo(outputVector->GetInformationObject(0));
 428       vtkDataSet *output(vtkDataSet::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
 429       std::vector<int> idsToKeep(this->Internal->getIdsToKeep());
 430       vtkDataSet *tryOnCell(FilterFamilies(input,idsToKeep,this->InsideOut));
 431       // first shrink the input
 432       output->ShallowCopy(tryOnCell);
 433       tryOnCell->Delete();
 434     }
 435   catch(INTERP_KERNEL::Exception& e)
 436     {
 437       std::cerr << "Exception has been thrown in vtkExtractCellType::RequestData : " << e.what() << std::endl;
 438       return 0;
 439     }
 440   return 1;
 441 }
 442
 443 int vtkExtractCellType::GetSILUpdateStamp()
 444 {
 445   return (int)this->SILTime;
 446 }
 447
 448 void vtkExtractCellType::PrintSelf(ostream& os, vtkIndent indent)
 449 {
 450   this->Superclass::PrintSelf(os, indent);
 451 }
 452
 453 int vtkExtractCellType::GetNumberOfGeoTypesArrays()
 454 {
 455   int ret(this->Internal->getNumberOfEntries());
 456   //std::cerr << "vtkExtractCellType::GetNumberOfGeoTypesArrays() -> " << ret << std::endl;
 457   return ret;
 458 }
 459
 460 const char *vtkExtractCellType::GetGeoTypesArrayName(int index)
 461 {
 462   const char *ret(this->Internal->getKeyOfEntry(index));
 463   //std::cerr << "vtkExtractCellType::GetGeoTypesArrayName(" << index << ") -> " << ret << std::endl;
 464   return ret;
 465 }
 466
 467 int vtkExtractCellType::GetGeoTypesArrayStatus(const char *name)
 468 {
 469   int ret((int)this->Internal->getStatusOfEntryStr(name));
 470   //std::cerr << "vtkExtractCellType::GetGeoTypesArrayStatus(" << name << ") -> " << ret << std::endl;
 471   return ret;
 472 }
 473
 474 void vtkExtractCellType::SetGeoTypesStatus(const char *name, int status)
 475 {
 476   //std::cerr << "vtkExtractCellType::SetGeoTypesStatus(" << name << "," << status << ")" << std::endl;
 477   if (GetNumberOfGeoTypesArrays()<1)
 478     return;
 479   this->Internal->setStatusOfEntryStr(name,(bool)status);
 480   if(std::string(name)==GetGeoTypesArrayName(GetNumberOfGeoTypesArrays()-1))
 481     {
 482       this->Modified();
 483       //this->Internal->printMySelf(std::cerr);
 484     }
 485 }