src/Plugins/SimpleMode/IO/vtkSimpleMode.cxx

   1 // Copyright (C) 2017  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 (EDF R&D)
  20
  21 #include "vtkSimpleMode.h"
  22
  23 #include "vtkAdjacentVertexIterator.h"
  24 #include "vtkDataSetSurfaceFilter.h"
  25 #include "vtkIntArray.h"
  26 #include "vtkCellData.h"
  27 #include "vtkPointData.h"
  28
  29 #include "vtkStreamingDemandDrivenPipeline.h"
  30 #include "vtkUnstructuredGrid.h"
  31 #include "vtkDataSet.h"
  32 #include  "vtkMultiBlockDataSet.h"
  33
  34 #include "vtkInformationStringKey.h"
  35 #include "vtkAlgorithmOutput.h"
  36 #include "vtkObjectFactory.h"
  37 #include "vtkMutableDirectedGraph.h"
  38 #include "vtkMultiBlockDataSet.h"
  39 #include "vtkDataSet.h"
  40 #include "vtkInformationVector.h"
  41 #include "vtkInformation.h"
  42 #include "vtkDataArraySelection.h"
  43 #include "vtkTimeStamp.h"
  44 #include "vtkInEdgeIterator.h"
  45 #include "vtkInformationDataObjectKey.h"
  46 #include "vtkExecutive.h"
  47 #include "vtkVariantArray.h"
  48 #include "vtkStringArray.h"
  49 #include "vtkDoubleArray.h"
  50 #include "vtkCharArray.h"
  51 #include "vtkUnsignedCharArray.h"
  52 #include "vtkDataSetAttributes.h"
  53 #include "vtkDemandDrivenPipeline.h"
  54 #include "vtkDataObjectTreeIterator.h"
  55 #include "vtkWarpScalar.h"
  56 #include "vtkWarpVector.h"
  57 #include "vtkMultiBlockDataGroupFilter.h"
  58 #include "vtkCompositeDataToUnstructuredGridFilter.h"
  59
  60 #ifdef WIN32
  61 #define _USE_MATH_DEFINES
  62 #include <math.h>
  63 #endif
  64
  65 #include <map>
  66 #include <deque>
  67 #include <sstream>
  68
  69 vtkStandardNewMacro(vtkSimpleMode);
  70
  71 static const char ZE_DISPLACEMENT_NAME1[]="@@ForReal?@@";
  72
  73 static const char ZE_DISPLACEMENT_NAME2[]="@@ForImag?@@";
  74
  75 static const char ZE_DISPLACEMENT_NAME3[]="MagnitudeOfCpxDisp";
  76
  77 static const double EPS=1e-12;
  78
  79 ///////////////////
  80
  81 class MZCException : public std::exception
  82 {
  83 public:
  84   MZCException(const std::string& s):_reason(s) { }
  85   virtual const char *what() const throw() { return _reason.c_str(); }
  86   virtual ~MZCException() throw() { }
  87 private:
  88   std::string _reason;
  89 };
  90
  91 vtkSmartPointer<vtkDoubleArray> ForceTo3Compo(vtkDoubleArray *arr)
  92 {
  93   if(!arr)
  94     return vtkSmartPointer<vtkDoubleArray>();
  95   int nbCompo(arr->GetNumberOfComponents()),nbTuples(arr->GetNumberOfTuples());
  96   if(nbCompo==3)
  97     {
  98       vtkSmartPointer<vtkDoubleArray> ret(arr);
  99       arr->Register(0);
 100       return ret;
 101     }
 102   if(nbCompo==6)
 103     {
 104       vtkSmartPointer<vtkDoubleArray> ret(vtkSmartPointer<vtkDoubleArray>::New());
 105       ret->SetNumberOfComponents(3);
 106       ret->SetNumberOfTuples(nbTuples);
 107       const double *srcPt(arr->Begin());
 108       double *destPt(ret->Begin());
 109       for(int i=0;i<nbTuples;i++,destPt+=3,srcPt+=6)
 110         std::copy(srcPt,srcPt+3,destPt);
 111       return ret;
 112     }
 113   throw MZCException("ForceTo3Compo : internal error ! 6 or 3 compo arrays expected !");
 114 }
 115
 116 std::vector< std::string > GetPossibleArrayNames(vtkDataSet *dataset)
 117 {
 118   if(!dataset)
 119     throw MZCException("The input dataset is null !");
 120   std::vector< std::string > ret;
 121   vtkPointData *att(dataset->GetPointData());
 122   for(int i=0;i<att->GetNumberOfArrays();i++)
 123     {
 124       vtkDataArray *locArr(att->GetArray(i));
 125       int nbComp(locArr->GetNumberOfComponents());
 126       if(nbComp!=3 && nbComp!=6)
 127         continue;
 128       std::string s(locArr->GetName());
 129       ret.push_back(s);
 130     }
 131   return ret;
 132 }
 133
 134 std::string FindTheBest(const std::vector<std::string>& arrNames, const std::string& key0, const std::string& key1)
 135 {
 136   std::string ret;
 137   char points(0);
 138   if(arrNames.empty())
 139     return ret;
 140   for(std::vector<std::string>::const_iterator it=arrNames.begin();it!=arrNames.end();it++)
 141     {
 142       char curNbPts(1);
 143       if((*it).find(key0,0)!=std::string::npos)
 144         curNbPts++;
 145       if((*it).find(key1,0)!=std::string::npos)
 146         curNbPts++;
 147       if(curNbPts>points)
 148         {
 149           points=curNbPts;
 150           ret=*it;
 151         }
 152     }
 153   return ret;
 154 }
 155
 156 std::string FindBestRealAmong(const std::vector<std::string>& arrNames)
 157 {
 158   static const char KEY1[]="DEPL";
 159   static const char KEY2[]="REEL";
 160   return FindTheBest(arrNames,KEY1,KEY2);
 161 }
 162
 163 std::string FindBestImagAmong(const std::vector<std::string>& arrNames)
 164 {
 165   static const char KEY1[]="DEPL";
 166   static const char KEY2[]="IMAG";
 167   return FindTheBest(arrNames,KEY1,KEY2);
 168 }
 169
 170 vtkUnstructuredGrid *ExtractInfo1(vtkInformationVector *inputVector)
 171 {
 172   vtkInformation *inputInfo(inputVector->GetInformationObject(0));
 173   vtkDataSet *input(0);
 174   vtkDataSet *input0(vtkDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
 175   vtkMultiBlockDataSet *input1(vtkMultiBlockDataSet::SafeDownCast(inputInfo->Get(vtkDataObject::DATA_OBJECT())));
 176   if(input0)
 177     input=input0;
 178   else
 179     {
 180       if(!input1)
 181         throw MZCException("Input dataSet must be a DataSet or single elt multi block dataset expected !");
 182       if(input1->GetNumberOfBlocks()!=1)
 183         throw MZCException("Input dataSet is a multiblock dataset with not exactly one block ! Use MergeBlocks or ExtractBlocks filter before calling this filter !");
 184       vtkDataObject *input2(input1->GetBlock(0));
 185       if(!input2)
 186         throw MZCException("Input dataSet is a multiblock dataset with exactly one block but this single element is NULL !");
 187       vtkDataSet *input2c(vtkDataSet::SafeDownCast(input2));
 188       if(!input2c)
 189         throw MZCException("Input dataSet is a multiblock dataset with exactly one block but this single element is not a dataset ! Use MergeBlocks or ExtractBlocks filter before calling this filter !");
 190       input=input2c;
 191     }
 192   if(!input)
 193     throw MZCException("Input data set is NULL !");
 194   vtkUnstructuredGrid *usgIn(vtkUnstructuredGrid::SafeDownCast(input));
 195   if(!usgIn)
 196     throw MZCException("Input data set is not an unstructured mesh ! This filter works only on unstructured meshes !");
 197   return usgIn;
 198 }
 199
 200 void ExtractInfo3(vtkDataSet *ds, const std::string& arrName, vtkDataArray *& arr, int& idx)
 201 {
 202   vtkPointData *att(ds->GetPointData());
 203   if(!att)
 204     throw MZCException("Input dataset has no point data attribute ! Impossible to move mesh !");
 205   vtkDataArray *zeArr(0);
 206   int i(0);
 207   for(;i<att->GetNumberOfArrays();i++)
 208     {
 209       vtkDataArray *locArr(att->GetArray(i));
 210       std::string s(locArr->GetName());
 211       if(s==arrName)
 212         {
 213           zeArr=locArr;
 214           break;
 215         }
 216     }
 217   if(!zeArr)
 218     {
 219       std::ostringstream oss;
 220       oss << "Impossible to locate the array called \"" << arrName << "\" used to move mesh !";
 221       throw MZCException(oss.str());
 222     }
 223   arr=zeArr;
 224   idx=i;
 225 }
 226
 227
 228 void ExtractInfo2(vtkDataSet *ds, const std::string& arrName, vtkDoubleArray *& arr)
 229 {
 230   vtkDataArray *zeArr(0);
 231   int dummy;
 232   ExtractInfo3(ds,arrName,zeArr,dummy);
 233   arr=vtkDoubleArray::SafeDownCast(zeArr);
 234   if(!arr)
 235     {
 236       std::ostringstream oss;
 237       oss << "Array called \"" << arrName << "\" has been located but this is NOT a float64 array !";
 238       throw MZCException(oss.str());
 239     }
 240   if(arr->GetNumberOfComponents()!=3 && arr->GetNumberOfComponents()!=6)
 241     {
 242       std::ostringstream oss;
 243       oss << "Float64 array called \"" << arrName << "\" has been located but this array has not exactly 3 or 6 components as it should !";
 244       throw MZCException(oss.str());
 245     }
 246   if(arr->GetNumberOfTuples()!=ds->GetNumberOfPoints())
 247     {
 248       std::ostringstream oss;
 249       oss << "Float64-1 components array called \"" << arrName << "\" has been located but the number of tuples is invalid ! Should be " << ds->GetNumberOfPoints() << " instead of " << arr->GetNumberOfTuples() << " !";
 250       throw MZCException(oss.str());
 251     }
 252 }
 253
 254 ////////////////////
 255
 256 class vtkSimpleMode::vtkSimpleModeInternal
 257 {
 258 public:
 259   vtkSimpleModeInternal():_surface(0) { }
 260   vtkPolyData *performConnection(vtkDataSet *ds);
 261   void setFieldForReal(const std::string& st) { _real=st; }
 262   std::string getFieldForReal() const { return _real; }
 263   ~vtkSimpleModeInternal();
 264 private:
 265   vtkDataSetSurfaceFilter *_surface;
 266 private:
 267   std::string _real;
 268 };
 269
 270 vtkPolyData *vtkSimpleMode::vtkSimpleModeInternal::performConnection(vtkDataSet *ds)
 271 {
 272   if(!_surface)
 273     {
 274       _surface=vtkDataSetSurfaceFilter::New();
 275       _surface->SetInputData(ds);
 276     }
 277   _surface->Update();
 278   return _surface->GetOutput();
 279 }
 280
 281 vtkSimpleMode::vtkSimpleModeInternal::~vtkSimpleModeInternal()
 282 {
 283   if(_surface)
 284     _surface->Delete();
 285 }
 286
 287 vtkSimpleMode::vtkSimpleMode():Factor(1.),AnimationTime(0.),Internal(new vtkSimpleMode::vtkSimpleModeInternal)
 288 {
 289   //this->SetInputArrayToProcess(0,0,0,vtkDataObject::FIELD_ASSOCIATION_POINTS,vtkDataSetAttributes::VECTORS);
 290 }
 291
 292 vtkSimpleMode::~vtkSimpleMode()
 293 {
 294   delete this->Internal;
 295 }
 296
 297 void vtkSimpleMode::SetInputArrayToProcess(int idx, int port, int connection, int ff, const char *name)
 298 {
 299   if(idx==0)
 300     this->Internal->setFieldForReal(name);
 301   vtkDataSetAlgorithm::SetInputArrayToProcess(idx,port,connection,ff,name);
 302 }
 303
 304 double GetOptimalRatioFrom(vtkUnstructuredGrid *dataset, vtkDoubleArray *array)
 305 {
 306   if(!dataset || !array)
 307     throw MZCException("The input dataset and or array is null !");
 308   vtkDataArray *coords(dataset->GetPoints()->GetData());
 309   vtkDoubleArray *coords2(vtkDoubleArray::SafeDownCast(coords));
 310   if(!coords2)
 311     throw MZCException("Input coordinates are not float64 !");
 312   int nbCompo(array->GetNumberOfComponents());
 313   if(coords2->GetNumberOfComponents()!=3 || (nbCompo!=3 && nbCompo!=6))
 314     throw MZCException("Input coordinates do not have 3 components as it should !");
 315   int nbPts(dataset->GetNumberOfPoints());
 316   const double *srcPt1(array->Begin());
 317   dataset->ComputeBounds();
 318   double *minmax1(dataset->GetBounds());
 319   double minmax2[3]={0.,0.,0.};
 320   for(int i=0;i<nbPts;i++,srcPt1+=nbCompo)
 321     {
 322       minmax2[0]=std::max(fabs(srcPt1[0]),minmax2[0]);
 323       minmax2[1]=std::max(fabs(srcPt1[1]),minmax2[1]);
 324       minmax2[2]=std::max(fabs(srcPt1[2]),minmax2[2]);
 325     }
 326   double maxDispDelta(*std::max_element(minmax2,minmax2+3));
 327   if(maxDispDelta<EPS)
 328     maxDispDelta=1.;
 329   for(int i=0;i<3;i++)
 330     minmax2[i]=minmax1[2*i+1]-minmax1[2*i];
 331   double maxGeoDelta(*std::max_element(minmax2,minmax2+3));
 332   if(maxDispDelta<EPS)
 333     maxDispDelta=1.;
 334   return maxGeoDelta/maxDispDelta;
 335 }
 336
 337 int vtkSimpleMode::FillOutputPortInformation( int vtkNotUsed(port), vtkInformation* info)
 338 {
 339   info->Set(vtkDataObject::DATA_TYPE_NAME(), "vtkPolyData");
 340   return 1;
 341 }
 342
 343 int vtkSimpleMode::RequestInformation(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 344 {
 345   //std::cerr << "########################################## vtkSimpleMode::RequestInformation ##########################################" << std::endl;
 346   try
 347     {
 348       if(this->Internal->getFieldForReal().empty())
 349         return 1;
 350       /*vtkUnstructuredGrid *usgIn(0);
 351       vtkDoubleArray *arr(0);
 352       ExtractInfo(inputVector[0],usgIn,this->Internal->getFieldForReal(),arr);
 353       std::vector<std::string> candidatesArrName(GetPossibleArrayNames(usgIn));
 354       //
 355       double ratio(GetOptimalRatioFrom(usgIn,arr));
 356       std::string optArrNameForReal(FindBestRealAmong(candidatesArrName));
 357       std::string optArrNameForImag(FindBestImagAmong(candidatesArrName));*/
 358       //std::cerr << ratio << std::endl;
 359       //std::cerr << optArrNameForReal << " * " << optArrNameForImag << std::endl;
 360     }
 361   catch(MZCException& e)
 362     {
 363       std::ostringstream oss;
 364       oss << "Exception has been thrown in vtkSimpleMode::RequestInformation : " << e.what() << std::endl;
 365       if(this->HasObserver("ErrorEvent") )
 366         this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
 367       else
 368         vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
 369       vtkObject::BreakOnError();
 370       return 0;
 371     }
 372   return 1;
 373 }
 374
 375 int vtkSimpleMode::RequestData(vtkInformation *request, vtkInformationVector **inputVector, vtkInformationVector *outputVector)
 376 {
 377   //std::cerr << "########################################## vtkSimpleMode::RequestData        ##########################################" << std::endl;
 378   try
 379     {
 380       vtkUnstructuredGrid *usgIn(0);
 381       usgIn=ExtractInfo1(inputVector[0]);
 382       //
 383       int nbPts(usgIn->GetNumberOfPoints());
 384       vtkPolyData *outSurface(this->Internal->performConnection(usgIn));
 385       vtkDoubleArray *arrRealBase(0);
 386       ExtractInfo2(outSurface,this->Internal->getFieldForReal(),arrRealBase);
 387       vtkSmartPointer<vtkDoubleArray> arrReal(ForceTo3Compo(arrRealBase));
 388       vtkSmartPointer<vtkDoubleArray> arr1(vtkSmartPointer<vtkDoubleArray>::New());
 389       arr1->SetName(ZE_DISPLACEMENT_NAME1);
 390       arr1->SetNumberOfComponents(3);
 391       arr1->SetNumberOfTuples(nbPts);
 392       double *ptToFeed1(arr1->Begin());
 393       const double *srcPt1(arrReal->Begin());
 394       double cst1(Factor*sin(AnimationTime*2*M_PI));
 395       std::transform(srcPt1,srcPt1+3*nbPts,ptToFeed1,std::bind2nd(std::multiplies<double>(),cst1));
 396       int idx1(outSurface->GetPointData()->AddArray(arr1));
 397       outSurface->GetPointData()->SetActiveAttribute(idx1,vtkDataSetAttributes::VECTORS);
 398       //
 399       //
 400       vtkSmartPointer<vtkWarpVector> ws(vtkSmartPointer<vtkWarpVector>::New());//vtkNew
 401       ws->SetInputData(outSurface);
 402       ws->SetScaleFactor(1.);
 403       ws->SetInputArrayToProcess(idx1,0,0,"vtkDataObject::FIELD_ASSOCIATION_POINTS",ZE_DISPLACEMENT_NAME1);
 404       ws->Update();
 405       vtkSmartPointer<vtkDataSet> ds(ws->GetOutput());
 406       ds->GetPointData()->RemoveArray(idx1);
 407       //
 408       int idx2(0);
 409       {
 410         vtkDataArray *dummy(0);
 411         ExtractInfo3(ds,this->Internal->getFieldForReal(),dummy,idx2);
 412       }
 413       ds->GetPointData()->SetActiveAttribute(idx2,vtkDataSetAttributes::SCALARS);
 414       //
 415       vtkInformation *outInfo(outputVector->GetInformationObject(0));
 416       vtkPolyData *output(vtkPolyData::SafeDownCast(outInfo->Get(vtkDataObject::DATA_OBJECT())));
 417       output->ShallowCopy(ds);
 418
 419     }
 420   catch(MZCException& e)
 421     {
 422       std::ostringstream oss;
 423       oss << "Exception has been thrown in vtkSimpleMode::RequestInformation : " << e.what() << std::endl;
 424       if(this->HasObserver("ErrorEvent") )
 425         this->InvokeEvent("ErrorEvent",const_cast<char *>(oss.str().c_str()));
 426       else
 427         vtkOutputWindowDisplayErrorText(const_cast<char *>(oss.str().c_str()));
 428       vtkObject::BreakOnError();
 429       return 0;
 430     }
 431   return 1;
 432 }
 433
 434 void vtkSimpleMode::PrintSelf(ostream& os, vtkIndent indent)
 435 {
 436   this->Superclass::PrintSelf(os, indent);
 437 }