src/libengrid/egvtkobject.h

   1 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   2 // +                                                                      +
   3 // + This file is part of enGrid.                                         +
   4 // +                                                                      +
   5 // + Copyright 2008-2014 enGits GmbH                                      +
   6 // +                                                                      +
   7 // + enGrid is free software: you can redistribute it and/or modify       +
   8 // + it under the terms of the GNU General Public License as published by +
   9 // + the Free Software Foundation, either version 3 of the License, or    +
  10 // + (at your option) any later version.                                  +
  11 // +                                                                      +
  12 // + enGrid is distributed in the hope that it will be useful,            +
  13 // + but WITHOUT ANY WARRANTY; without even the implied warranty of       +
  14 // + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the        +
  15 // + GNU General Public License for more details.                         +
  16 // +                                                                      +
  17 // + You should have received a copy of the GNU General Public License    +
  18 // + along with enGrid. If not, see <http://www.gnu.org/licenses/>.       +
  19 // +                                                                      +
  20 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  21 #ifndef egvtkobject_H
  22 #define egvtkobject_H
  23
  24 class EgVtkObject;
  25
  26 class BezierTriangle;
  27
  28 #include "engrid.h"
  29 #include "utilities.h"
  30 #include "boundarycondition.h"
  31
  32 #include <vtkUnstructuredGrid.h>
  33 #include <vtkPolyData.h>
  34 #include <vtkPointData.h>
  35 #include <vtkCellData.h>
  36 #include <vtkLongArray.h>
  37 #include <vtkDoubleArray.h>
  38 #include <vtkXMLUnstructuredGridWriter.h>
  39 #include <vtkCellLocator.h>
  40 #include <vtkIdList.h>
  41
  42 #include <QSettings>
  43 #include <QSet>
  44 #include <QVector>
  45
  46 #define EG_SIMPLE_VERTEX         0
  47 #define EG_FEATURE_EDGE_VERTEX   1
  48 #define EG_BOUNDARY_EDGE_VERTEX  2
  49 #define EG_FEATURE_CORNER_VERTEX 3
  50 #define EG_FIXED_VERTEX          4
  51
  52 class EgVtkObject
  53 {
  54
  55 public: // data-types
  56
  57   typedef const QVector<vtkIdType>&     l2g_t;
  58   typedef const QVector<int>&           g2l_t;
  59   typedef const QVector<QVector<int> >& l2l_t;
  60
  61 private: // methods
  62
  63   void addToC2C
  64     (
  65       vtkIdType               id_cell,
  66       QVector<int>           &_cells,
  67       QVector<QVector<int> > &c2c,
  68       int                     j,
  69       vtkIdList              *nds,
  70       vtkIdList              *cls,
  71       vtkUnstructuredGrid    *grid
  72     );
  73
  74   void addToN2N
  75     (
  76       QVector<QSet<int> > &n2n,
  77       int                  n1,
  78       int                  n2
  79     );
  80
  81   void createNodeField(vtkUnstructuredGrid *grid, QString field_name, QString type_name, int Nnodes, bool overwrite = false);
  82   void createCellField(vtkUnstructuredGrid *grid, QString field_name, QString type_name, int Ncells, bool overwrite = false);
  83
  84   QString getXmlSection(QString name);
  85
  86 protected: // attributes
  87
  88   QSet<int> m_BoundaryCodes;
  89   static int DebugLevel;
  90
  91 protected: // methods
  92
  93   /**
  94    * if key=value pair not found in settings file, write it + read key value from settings file and assign it to variable
  95    * Version for int variables
  96    */
  97   int getSet(QString group, QString key, int value, int& variable);
  98
  99   /**
 100    * if key=value pair not found in settings file, write it + read key value from settings file and assign it to variable
 101    * Version for double variables
 102    */
 103   double getSet(QString group, QString key, double value, double& variable);
 104
 105   /**
 106    * if key=value pair not found in settings file, write it + read key value from settings file and assign it to variable
 107    * Version for bool variables
 108    */
 109   bool getSet(QString group, QString key, bool value, bool& variable);
 110
 111   /**
 112    * if key=value pair not found in settings file, write it + read key value from settings file and assign it to variable
 113    * Version for string variables
 114    */
 115   QString getSet(QString group, QString key, QString value, QString& variable);
 116
 117   /**
 118    * if key=value pair not found in settings file, write it + read key value from settings file and assign it to variable
 119    * Version for string variables.
 120    */
 121   QString getSet(QString group, QString key, QString value, QString& variable, int type);
 122
 123   template <typename T>
 124   bool getXmlSetting(QString key, QString xml_section, T& value);
 125
 126   /**
 127    * Update the cell index array.
 128    */
 129   void UpdateCellIndex(vtkUnstructuredGrid *grid);
 130
 131   /**
 132    * Update the point index array.
 133    */
 134   void UpdateNodeIndex(vtkUnstructuredGrid *grid);
 135
 136   /**
 137    * Compute normal vectors on nodes and cells of a subset of a grid.
 138    * The paramters nodes and cells must be consistent; this means the nodes
 139    * represent exactly (not more, not less) the nodes forming the cells.
 140    * @param cell_normals On return, this will contain the cell normals (same order as cells)
 141    * @param node_normals On return, this will contain the cell normals (same order as cells)
 142    * @param cells        The cells to compute the normals of
 143    * @param nodes        The nodes to compute the normals of
 144    * @param grid         The grid to operate on
 145    */
 146   void computeNormals(QVector<vec3_t> &cell_normals, QVector<vec3_t> &node_normals, QVector<vtkIdType> &cells, QVector<vtkIdType>  &nodes, vtkUnstructuredGrid *grid);
 147
 148   /**
 149    * Create a mapping from global node indices to the indeces of a subset of nodes.
 150    * @param nodes  The subset of nodes.
 151    * @param _nodes On return, this will contain the mapping.
 152    * @param grid   The grid to operate on.
 153    */
 154   void createNodeMapping(QVector<vtkIdType> &nodes, QVector<int> &_nodes, vtkUnstructuredGrid *grid);
 155
 156   /**
 157    * Create a mapping from global cell indices to the indices of a subset of cells.
 158    * @param cells  The subset of cells.
 159    * @param _cells On return, this will contain the mapping.
 160    * @param grid   The grid to operate on.
 161    */
 162   void createCellMapping(QVector<vtkIdType> &cells, QVector<int> &_cells, vtkUnstructuredGrid *grid);
 163
 164   /**
 165    * Create a node to boundary condition ("cell_code") mapping.
 166    * Only non-zero boundary conditions will be considered.
 167    * @param bcs  On return, this will hold the codes of all boundary elements that are
 168    *             attached to a node.
 169    * @param grid The grid to operate on.
 170    */
 171   void createNodeToBcMapping(QVector<QSet<int> > &bcs, vtkUnstructuredGrid *grid);
 172
 173   /**
 174    * Create a node to cell structure for a given set of cells and nodes.
 175    * This creates a vector of sets which might have performance issues.
 176    * @param cells  the subset of cells
 177    * @param nodes  the subset of nodes
 178    * @param _nodes the reverse mapping for the nodes
 179    * @param n2c    On return, this will hold the node to cell structure
 180    * @param grid   The grid to operate on
 181    */
 182   void createNodeToCell(QVector<vtkIdType> &cells, QVector<vtkIdType> &nodes, QVector<int> &_nodes, QVector<QSet<int> >  &n2c, vtkUnstructuredGrid *grid);
 183
 184   /**
 185    * Create a node to cell structure for a given set of cells and nodes.
 186    * This creates a vector of vectors.
 187    * @param cells  the subset of cells
 188    * @param nodes  the subset of nodes
 189    * @param _nodes the reverse mapping for the nodes
 190    * @param n2c    On return, this will hold the node to cell structure
 191    * @param grid   The grid to operate on
 192    */
 193   void createNodeToCell(QVector<vtkIdType> &cells, QVector<vtkIdType> &nodes, QVector<int> &_nodes, QVector<QVector<int> > &n2c, vtkUnstructuredGrid *grid);
 194
 195   /**
 196    * Create a node to node structure for a given set of cells and nodes.
 197    * This creates a vector of sets which might have performance issues.
 198    * @param cells  the subset of cells
 199    * @param nodes  the subset of nodes
 200    * @param _nodes the reverse mapping for the nodes
 201    * @param n2n    On return, this will hold the node to node structure
 202    * @param grid   The grid to operate on
 203    */
 204   void createNodeToNode(QVector<vtkIdType> &cells, QVector<vtkIdType> &nodes, QVector<int> &_nodes, QVector<QSet<int> > &n2n, vtkUnstructuredGrid *grid);
 205
 206   /**
 207    * Create a node to node structure for a given set of cells and nodes.
 208    * This creates a vector of vectors.
 209    * @param cells  the subset of cells
 210    * @param nodes  the subset of nodes
 211    * @param _nodes the reverse mapping for the nodes
 212    * @param n2n    On return, this will hold the node to node structure
 213    * @param grid   The grid to operate on
 214    */
 215   void createNodeToNode(QVector<vtkIdType> &cells, QVector<vtkIdType> &nodes, QVector<int> &_nodes, QVector<QVector<int> > &n2n, vtkUnstructuredGrid *grid);
 216
 217   /**
 218    * Extract the nodes which are part of a given set of cells.
 219    * @param cells the subset of cells
 220    * @param nodes On return, this will contain the nodes that correspond to the subset of cells
 221    * @param grid  The grid to operate on
 222    */
 223   template <class C>
 224   void getNodesFromCells(const C &cells, QVector<vtkIdType> &nodes, vtkUnstructuredGrid *grid);
 225
 226   /**
 227    * Check if a cell is a volume cell.
 228    * @param cellId The id fof the cell in question
 229    * @param grid   The grid to operate on
 230    * @return true if the cell represents a volume and false if not
 231    */
 232   bool isVolume(vtkIdType id_cell, vtkUnstructuredGrid *grid);
 233
 234
 235   /**
 236    * Check if a cell is a surface cell.
 237    * @param cellId The id fof the cell in question
 238    * @param grid   The grid to operate on
 239    * @return true if the cell represents a surface and false if not
 240    */
 241   bool isSurface(vtkIdType id_cell, vtkUnstructuredGrid *grid);
 242
 243   /**
 244    * Get all volume cells of a grid.
 245    * @param cells On return this will hold the Ids of all volume cells.
 246    * @param grid  The grid to operate on.
 247    */
 248   void getAllVolumeCells(QVector<vtkIdType> &cells, vtkUnstructuredGrid *grid);
 249
 250   /**
 251    * Get all cells of a grid.
 252    * @param cells On return this will hold the Ids of all cells.
 253    * @param grid  The grid to operate on.
 254    */
 255   void getAllCells(QVector<vtkIdType> &cells, vtkUnstructuredGrid *grid);
 256
 257   /**
 258    * Get all cells of a grid and a specific type.
 259    * @param type  The type of the cells (e.g. VTK_TETRA, VTK_TRIANGLE, etc.)
 260    * @param cells On return this will hold the Ids of all cells.
 261    * @param grid  The grid to operate on.
 262    */
 263   void getAllCellsOfType(vtkIdType type, QVector<vtkIdType> &cells, vtkUnstructuredGrid *grid);
 264
 265   /**
 266    * Get all surface cells of a grid.
 267    * @param cells On return this will hold the Ids of all surface cells.
 268    * @param grid  The grid to operate on.
 269    */
 270   void getAllSurfaceCells(QVector<vtkIdType>  &cells, vtkUnstructuredGrid *grid);
 271
 272   /**
 273    * Get all surface cells of a grid with a specific boundary condition.
 274    * @param bcs   The set of boundary conditions
 275    * @param cells On return this will hold the Ids of the surface cells.
 276    * @param grid  The grid to operate on.
 277    */
 278   template <typename C>
 279   void getSurfaceCells(const C &bcs, QVector<vtkIdType> &cells, vtkUnstructuredGrid *grid);
 280
 281   /**
 282    * Create a cell neighbourship list for a subset grid.
 283    * This has been implemented using VTK's vtkCellLinks structures.
 284    * @param cells the subset of cells
 285    * @param c2c   On return this will hold the neighbourship list
 286    * @param grid  The grid to operate on.
 287    */
 288   void createCellToCell(QVector<vtkIdType> &cells, QVector<QVector<int> > &c2c, vtkUnstructuredGrid *grid);
 289
 290   /**
 291    * Insert a subset of a grid into a vtkPolyData structure.
 292    * This is can be used in order to make use of many readily available
 293    * operations within VTK; one example is smoothing.
 294    * Cell index and node index arrays will be created and passed to the
 295    * poly data structure. Thus any purely geometric effect (no topology change)
 296    * can be directly reintroduced into the vtkUnstructuredGrid.
 297    * @param cells the subset of cells
 298    * @param pdata the vtkPolyData to add the nodes and cells to
 299    * @param grid  The grid to operate on.
 300    */
 301   void addToPolyData(QVector<vtkIdType>  &cells, vtkPolyData *pdata, vtkUnstructuredGrid *grid);
 302
 303   /**
 304    * Copy the attributes from an input to an output cell.
 305    * @param old_grid the input grid
 306    * @param oldId    the existing input cell
 307    * @param new_grid the output grid
 308    * @param newId    the new output cell
 309    */
 310   void copyCellData(vtkUnstructuredGrid *old_grid, vtkIdType oldId, vtkUnstructuredGrid *new_grid, vtkIdType newId);
 311
 312   /**
 313    * Copy the attributes from an input to an output node.
 314    * @param old_grid the input grid
 315    * @param oldId    the existing input node
 316    * @param new_grid the output grid
 317    * @param newId    the new output node
 318    */
 319   void copyNodeData(vtkUnstructuredGrid *old_grid, vtkIdType oldId, vtkUnstructuredGrid *new_grid, vtkIdType newId);
 320
 321   /**
 322    * Create the basic fields on a given grid.
 323    * Care should be taken with the overwrite parameter; if it is set to <i>false</i>
 324    * and the field does not have the correct size it can lead to <i>segmentation faults</i>.
 325    * @param Ncells the number of output cells
 326    * @param Nnodes the number of output nodes
 327    * @param overwrite f set to true existing fields will be re-created
 328    */
 329   void createBasicFields(vtkUnstructuredGrid *grid, vtkIdType num_cells = -1, vtkIdType num_nodes = -1, bool overwrite = false);
 330
 331   /**
 332    * Create the basic cell fields on a given grid.
 333    * Care should be taken with the overwrite parameter; if it is set to <i>false</i>
 334    * and the field does not have the correct size it can lead to <i>segmentation faults</i>.
 335    * @param Ncells the number of output cells
 336    * @param overwrite f set to true existing fields will be re-created
 337    */
 338   void createBasicCellFields(vtkUnstructuredGrid *grid, vtkIdType Ncells, bool overwrite = false);
 339
 340   /**
 341    * Create the basic node fields on a given grid.
 342    * Care should be taken with the overwrite parameter; if it is set to <i>false</i>
 343    * and the field does not have the correct size it can lead to <i>segmentation faults</i>.
 344    * @param Nnodes the number of output nodes
 345    * @param overwrite f set to true existing fields will be re-created
 346    */
 347   void createBasicNodeFields(vtkUnstructuredGrid *grid, vtkIdType Nnodes, bool overwrite = false);
 348
 349   /**
 350    * Allocate memory for a grid. This method will also create the basic
 351    * attribute fields (e.g. "cell_code").
 352    * @param grid   the grid for which to allocate memory
 353    * @param Ncells the number of output cells
 354    * @param Nnodes the number of output nodes
 355    * @param create_fields flag to determine if node and cell data shall be created
 356    */
 357   void allocateGrid(vtkUnstructuredGrid *grid, vtkIdType Ncells, vtkIdType Nnodes, bool create_fields = true);
 358
 359   /**
 360    * Get the names of all node (point) attributes (fields) of a VTK grid
 361    * @param field_names On return this vector will contain the names of all fields
 362    * @param grid the grid
 363    */
 364   void getAllNodeDataNames(QVector<QString> &field_names, vtkUnstructuredGrid *grid);
 365
 366   /**
 367    * Get the names of all cell attributes (fields) of a VTK grid
 368    * @param field_names On return this vector will contain the names of all fields
 369    * @param grid the grid
 370    */
 371   void getAllCellDataNames(QVector<QString> &field_names, vtkUnstructuredGrid *grid);
 372
 373   /**
 374    * Compute the intersection of two Q containers.
 375    * This will return a set.
 376    * @param set1 the first container
 377    * @param set2 the second container
 378    * @param inters on return this will hold the intersection
 379    */
 380   template <class C1, class C2>
 381   void qcontIntersection(const C1& c1, const C2& c2, QSet<typename C1::value_type> &inters);
 382
 383   /**
 384    * Compute the intersection of two Q containers.
 385    * This will return a vector.
 386    * @param set1 the first container
 387    * @param set2 the second container
 388    * @param inters on return this will hold the intersection
 389    */
 390   template <class C1, class C2>
 391   void qcontIntersection(const C1& c1, const C2& c2, QVector<typename C1::value_type> &inters);
 392
 393   /**
 394    * Compute the centre of a cell
 395    * @param grid the grid to use
 396    * @param id_cell the cell of which to compute the centre
 397    * @return the centre of the cell
 398    */
 399   vec3_t cellCentre(vtkUnstructuredGrid *grid, vtkIdType id_cell);
 400
 401   /**
 402    * Compute the angle between two faces.
 403    * @param grid the grid to use
 404    * @param id_face1 index of the first face
 405    * @param id_face2 index of the second face
 406    * @return the angle between the faces (M_PI for a flat surface)
 407    */
 408   double faceAngle(vtkUnstructuredGrid* grid, vtkIdType id_face1, vtkIdType id_face2);
 409
 410   /**
 411    * Get the cells of a grid that are not part of a given set of cells.
 412    * @param grid the grid to use
 413    * @param cells the given set of cells
 414    * @param rest_cells on return this will hold the rest of all cells of the grid (not part of cells)
 415    */
 416   void getRestCells(vtkUnstructuredGrid *grid, const QVector<vtkIdType> &cells, QVector<vtkIdType> &rest_cells);
 417
 418   /**
 419    * Find the corresponding volume cell of a surface cell
 420    * @param grid the grid to use
 421    * @param id_surf the id of the surface cell
 422    * @param n2n the node to cell structure for this grid
 423    * @return the id of the corresponding volume cell (or -1 if not found)
 424    */
 425   vtkIdType findVolumeCell(vtkUnstructuredGrid *grid, vtkIdType id_surf, g2l_t _nodes, l2g_t cells, g2l_t _cells, l2l_t n2c);
 426
 427   /**
 428    * @brief Copy a cell from one grid to another.
 429    * @param src the source grid
 430    * @param id_cell the cell index within the source grid
 431    * @param dst the destination grid
 432    * @return the index of the cell in the destination grid
 433    */
 434   vtkIdType copyCell(vtkUnstructuredGrid* src, vtkIdType id_cell, vtkUnstructuredGrid* dst, vtkIdType offset = 0);
 435
 436   /**
 437    * @brief Copy a cell from one grid to another and translate node indices.
 438    * @param src the source grid
 439    * @param id_cell the cell index within the source grid
 440    * @param dst the destination grid
 441    * @param src2dst a generic container containing the mapping from src to destination index
 442    * @return the index of the cell in the destination grid
 443    */
 444   template <typename C>
 445   vtkIdType copyCell(vtkUnstructuredGrid* src, vtkIdType id_cell, vtkUnstructuredGrid* dst, const C& src2dst);
 446
 447   /**
 448    * Copy "src" grid to "dst" grid. Allocate "dst" so that it fits the data of "src".
 449    * @param src a pointer to the source grid
 450    * @param dst a pointer to the destination grid
 451    */
 452   void makeCopy(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst, bool copy_data = true);
 453
 454   /**
 455    * Copy a part of "src" grid to "dst" grid. Allocate "dst" so that it fits the data to be copied.
 456    * @param src a pointer to the source grid
 457    * @param dst a pointer to the destination grid
 458    * @param cells a container with the cells to be copied
 459    */
 460   template <class C>
 461   void makeCopy(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst, const C &cells);
 462
 463   /**
 464    * Copy "src" grid to "dst" grid. DO NOT allocate "dst" so that it fits the data of "src".
 465    * Allocation is left for the user to do.
 466    * @param src a pointer to the source grid
 467    * @param dst a pointer to the destination grid
 468    */
 469   void makeCopyNoAlloc(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst);
 470
 471   /**
 472    * Change the orientation of a face.
 473    * @param grid the grid to use
 474    * @param id_face the id of the face to change
 475    */
 476   void reorientateFace(vtkUnstructuredGrid *grid, vtkIdType id_face);
 477
 478   /**
 479    * Reset face orientation to original orientation.
 480    * @param grid the grid with the faces
 481    */
 482   void resetOrientation(vtkUnstructuredGrid *grid);
 483
 484   void createIndices(vtkUnstructuredGrid *grid);
 485
 486   /**
 487    * Get the boundary condition of a boundary code.
 488    * @param bc the boundary code
 489    * @return the boundary condition
 490    */
 491   BoundaryCondition getBC(int bc);
 492
 493   /**
 494    * Save the subgrid defined by cls from grid.
 495    * @param grid The source grid
 496    * @param cls The cells to extract
 497    * @param file_name The file to save to
 498    */
 499   template <class C>
 500   void writeCells(vtkUnstructuredGrid *grid, const C &cls, QString file_name);
 501
 502   /**
 503    * Get the SubGrid defined by cls from grid. The function takes care of allocation for SubGrid.
 504    * @param grid The source grid
 505    * @param cls The cells to extract
 506    * @param SubGrid The SubGrid to create
 507    */
 508   template <class C>
 509   void getSubGrid(vtkUnstructuredGrid *grid, const C &cls, vtkUnstructuredGrid *SubGrid);
 510
 511
 512   /**
 513    * Save grid to file filename.
 514    * @param grid The source grid
 515    * @param filename Name of the file to save to
 516    */
 517   void writeGrid(vtkUnstructuredGrid *grid, QString filename);
 518
 519   /**
 520    * Get a file name without extension.
 521    * @param file_name the full name (with extension)
 522    * @return the name without the extension
 523    */
 524   QString stripFromExtension(QString file_name);
 525
 526   /**
 527    * Get the extension of a file name
 528    * @param file_name the full name (with extension)
 529    * @return the extension
 530    */
 531   QString getExtension(QString file_name);
 532
 533   /**
 534    * Get a face of a cell
 535    * @param grid the unstructured grid
 536    * @param id_cell the index of the cell
 537    * @param i_face the index of the face within the cell
 538    * @param ids on return this will contain the nodes of the face
 539    */
 540   void getFaceOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_face, QVector<vtkIdType> &ids);
 541
 542   /**
 543    * Get the normal of a face of a volume cell.
 544    * @param grid the unstructured grid
 545    * @param id_cell the index of the cell
 546    * @param i_face the index of the face within the cell
 547    * @return the normal vector (absolute value corresponds to the area)
 548    */
 549   vec3_t getNormalOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_face);
 550
 551   /**
 552    * Get the centre of a face of a volume cell.
 553    * @param grid the unstructured grid
 554    * @param id_cell the index of the cell
 555    * @param i_face the index of the face within the cell
 556    * @return the centre
 557    */
 558   vec3_t getCentreOfCellFace(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_face);
 559
 560   /**
 561    * Get an edge of a face/cell
 562    * @param grid the unstructured grid
 563    * @param id_cell the index of the cell
 564    * @param i_edge the index of the edge within the cell
 565    * @param ids on return this will contain the nodes of the edge
 566    */
 567   void getEdgeOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_edge, QVector<vtkIdType> &ids);
 568
 569   /**
 570    * Get all boundary codes fo a grid.
 571    * @param grid the grid to extract the boundaru codes from
 572    * @return a set with all boundary codes
 573    */
 574   QSet<int> getAllBoundaryCodes(vtkUnstructuredGrid *grid);
 575
 576   /**
 577    * Check if a cell (face) contains a given node
 578    * @param grid the unstructured grid
 579    * @param id_cell the id of the cell to investigate
 580    * @param id_node the id of the required node
 581    * @return true if id_cell contains id_node
 582    */
 583   bool cellContainsNode(vtkUnstructuredGrid *grid, vtkIdType id_cell, vtkIdType id_node);
 584
 585   /**
 586    * Get all node indices which are shared by two cells.
 587    * These cells can be surface or volume cells; also a combination
 588    * of a volume and a surface cell is possible.
 589    * @param grid the grid to use
 590    * @param id_cell1 index of the first cell
 591    * @param id_cell2 index of the second cell
 592    * @param cont a generic Qt container which will hold the shared node indices on return
 593    */
 594   template <typename C>
 595   void sharedNodesOfCells(vtkUnstructuredGrid* grid, vtkIdType id_cell1, vtkIdType id_cell2, C& cont);
 596
 597   /**
 598    * @brief Get all nodes of a cell from a vtkUnstructuredGrid
 599    * This methods collects all nodes of a cell in a generic Qt container.
 600    * It can be used uniformly for VTK_POLYHEDRON cells and standard cells.
 601    * @param grid the grid to use
 602    * @param id_cell the cell index
 603    * @param cont a generic Qt container which will hold node indices on return
 604    */
 605   template <typename C>
 606   void getPointsOfCell(vtkUnstructuredGrid* grid, vtkIdType id_cell, C& cont);
 607
 608   template <class C> void createPolyData(const C &x, vtkPolyData *poly_data, bool closed_loop = false);
 609   void createPolyDataC2C(vtkPolyData *poly_data, QVector<QVector<vtkIdType> > &c2c);
 610   void createPolyDataN2C(vtkPolyData *poly_data, QVector<QSet<vtkIdType> > &n2c);
 611   void createPolyDataN2N(vtkPolyData *poly_data, QVector<QSet<vtkIdType> > &n2n);
 612   template <class C> double convexRatio(const C &x, vec3_t n_plane, bool closed_loop = false);
 613
 614   template <class C> void invertQContainer(C &cont);
 615
 616
 617 public: // methods
 618
 619   EgVtkObject() { DebugLevel = 0; }
 620
 621   void setBoundaryCodes(const QSet<int> &bcs);
 622   QSet<int> getBoundaryCodes();
 623   void setDebugLevel(int a_DebugLevel) { DebugLevel = a_DebugLevel; }
 624
 625   bool saveGrid( vtkUnstructuredGrid* a_grid, QString file_name );
 626
 627   vtkIdType addGrid(vtkUnstructuredGrid *main_grid, vtkUnstructuredGrid *grid_to_add, vtkIdType offset);
 628
 629   template <typename C>
 630   void vtkIdList2QContainer(vtkIdList *id_list, C &cont);
 631
 632   void checkGridConsitency(vtkUnstructuredGrid *grid);
 633
 634   template <typename C1, typename C2>
 635   void triangulatePolygon(vtkUnstructuredGrid *grid, const C1 &polygon, C2 &triangles);
 636
 637 private:
 638
 639   void addVtkTypeInfo(vtkUnstructuredGrid* a_grid); ///< Add VTK type information to the grid (useful for visualisation with ParaView).
 640
 641 };
 642
 643 //End of class EgVtkObject
 644
 645 template <class C1, class C2>
 646 void EgVtkObject::qcontIntersection(const C1& c1, const C2& c2, QSet<typename C1::value_type> &inters)
 647 {
 648   inters.clear();
 649   foreach (typename C1::value_type t1, c1) {
 650     foreach (typename C2::value_type t2, c2) {
 651       if (t1 == t2) {
 652         inters.insert(t1);
 653       }
 654     }
 655   }
 656 }
 657
 658 template <class C1, class C2>
 659 void EgVtkObject::qcontIntersection(const C1& c1, const C2& c2, QVector<typename C1::value_type> &inters)
 660 {
 661   QSet<typename C1::value_type> inters_set;
 662   qcontIntersection(c1, c2, inters_set);
 663   inters.resize(inters_set.size());
 664   qCopy(inters_set.begin(), inters_set.end(), inters.begin());
 665 }
 666
 667 template <class C>
 668 void EgVtkObject::getSubGrid(vtkUnstructuredGrid *grid, const C &cls, vtkUnstructuredGrid *sub_grid)
 669 {
 670   createIndices(grid);
 671   QVector<vtkIdType> cells;
 672   QVector<vtkIdType> nodes;
 673   cells.resize(cls.size());
 674   qCopy(cls.begin(), cls.end(), cells.begin());
 675   getNodesFromCells(cells, nodes, grid);
 676   allocateGrid(sub_grid, cells.size(), nodes.size());
 677   vtkIdType id_new_node = 0;
 678   QVector<vtkIdType> old2new(grid->GetNumberOfPoints(), -1);
 679   foreach (vtkIdType id_node, nodes) {
 680     vec3_t x;
 681     grid->GetPoint(id_node, x.data());
 682     sub_grid->GetPoints()->SetPoint(id_new_node, x.data());
 683     old2new[id_node] = id_new_node;
 684     copyNodeData(grid, id_node, sub_grid, id_new_node);
 685     ++id_new_node;
 686   }
 687   foreach (vtkIdType id_cell, cells) {
 688     vtkIdType id_new_cell = copyCell(grid, id_cell, sub_grid, old2new);
 689     copyCellData(grid, id_cell, sub_grid, id_new_cell);
 690   }
 691 }
 692
 693 template <class C>
 694 void EgVtkObject::writeCells(vtkUnstructuredGrid *grid, const C &cls, QString file_name)
 695 {
 696   qDebug()<<"Saving cells from grid as "<<file_name;
 697
 698   EG_VTKSP(vtkUnstructuredGrid,SubGrid);
 699   getSubGrid(grid,cls,SubGrid);
 700
 701   EG_VTKSP(vtkXMLUnstructuredGridWriter,vtu);
 702   vtu->SetFileName(qPrintable(file_name));
 703   vtu->SetDataModeToBinary();
 704   vtu->SetInputData(SubGrid);
 705   vtu->Write();
 706 }
 707
 708 template <class C>
 709 void EgVtkObject::getNodesFromCells(const C& cells, QVector<vtkIdType>  &nodes, vtkUnstructuredGrid *grid)
 710 {
 711   QSet<vtkIdType> ex_nodes;
 712   vtkIdType id_cell;
 713   foreach(id_cell, cells) {
 714     QList<vtkIdType> pts;
 715     getPointsOfCell(grid, id_cell, pts);
 716     foreach (vtkIdType id_node, pts) {
 717       if (id_node >= grid->GetNumberOfPoints()) {
 718         EG_BUG;
 719       }
 720       ex_nodes.insert(id_node);
 721     }
 722   }
 723   nodes.resize(ex_nodes.size());
 724   {
 725     int j = 0;
 726     vtkIdType i;
 727     foreach(i,ex_nodes) {
 728       nodes[j] = i;
 729       ++j;
 730     }
 731   }
 732 }
 733
 734 inline vtkIdType EgVtkObject::copyCell(vtkUnstructuredGrid *src, vtkIdType id_cell, vtkUnstructuredGrid *dst, vtkIdType offset)
 735 {
 736   EG_VTKSP(vtkIdList, stream);
 737   vtkIdType type_cell = src->GetCellType(id_cell);
 738   vtkIdType id_new_cell = -1;
 739   if (type_cell == VTK_POLYHEDRON) {
 740     src->GetFaceStream(id_cell, stream);
 741     int id = 1;
 742     for (int i = 0; i < stream->GetId(0); ++i) {
 743       int num_pts = stream->GetId(id);
 744       ++id;
 745       for (int j = 0; j < num_pts; ++j) {
 746         stream->SetId(id, stream->GetId(id) + offset);
 747         ++id;
 748       }
 749     }
 750     id_new_cell = dst->InsertNextCell(type_cell, stream);
 751   } else {
 752     src->GetCellPoints(id_cell, stream);
 753     for (int i = 0; i < stream->GetNumberOfIds(); ++i) {
 754       stream->SetId(i, stream->GetId(i) + offset);
 755     }
 756     id_new_cell = dst->InsertNextCell(type_cell, stream);
 757   }
 758   return id_new_cell;
 759 }
 760
 761 template <typename C>
 762 vtkIdType EgVtkObject::copyCell(vtkUnstructuredGrid *src, vtkIdType id_cell, vtkUnstructuredGrid *dst, const C& src2dst)
 763 {
 764   EG_VTKSP(vtkIdList, stream);
 765   vtkIdType type_cell = src->GetCellType(id_cell);
 766   vtkIdType id_new_cell = -1;
 767   if (type_cell == VTK_POLYHEDRON) {
 768     src->GetFaceStream(id_cell, stream);
 769     int id = 1;
 770     for (int i = 0; i < stream->GetId(0); ++i) {
 771       int num_pts = stream->GetId(id);
 772       ++id;
 773       for (int j = 0; j < num_pts; ++j) {
 774         stream->SetId(id, src2dst[stream->GetId(id)]);
 775         ++id;
 776       }
 777     }
 778     id_new_cell = dst->InsertNextCell(type_cell, stream);
 779   } else {
 780     src->GetCellPoints(id_cell, stream);
 781     for (int i = 0; i < stream->GetNumberOfIds(); ++i) {
 782       if (src2dst[stream->GetId(i)] < 0) {
 783         EG_BUG;
 784       }
 785       if (src2dst[stream->GetId(i)] >= dst->GetNumberOfPoints()) {
 786         EG_BUG;
 787       }
 788       stream->SetId(i, src2dst[stream->GetId(i)]);
 789     }
 790     id_new_cell = dst->InsertNextCell(type_cell, stream);
 791   }
 792   return id_new_cell;
 793 }
 794
 795 template <class C>
 796 void EgVtkObject::makeCopy(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst, const C& cells)
 797 {
 798   QVector<vtkIdType> nodes;
 799   getNodesFromCells(cells, nodes, src);
 800   allocateGrid(dst, cells.size(), nodes.size());
 801   vtkIdType id_new_node = 0;
 802   QVector<vtkIdType> old2new(src->GetNumberOfPoints(), -1);
 803   foreach (vtkIdType id_node, nodes) {
 804     vec3_t x;
 805     src->GetPoints()->GetPoint(id_node, x.data());
 806     dst->GetPoints()->SetPoint(id_new_node, x.data());
 807     copyNodeData(src, id_node, dst, id_new_node);
 808     old2new[id_node] = id_new_node;
 809     ++id_new_node;
 810   }
 811   foreach (vtkIdType id_cell, cells) {
 812     vtkIdType id_new_cell = copyCell(src, id_cell, dst, old2new);
 813     copyCellData(src, id_cell, dst, id_new_cell);
 814   }
 815 }
 816
 817 template <class C>
 818 void EgVtkObject::createPolyData(const C &x, vtkPolyData *poly_data, bool closed_loop)
 819 {
 820   int N = x.size();
 821   if (closed_loop) {
 822     --N;
 823   }
 824   EG_VTKSP(vtkPoints, points);
 825   points->SetNumberOfPoints(N);
 826   QVector<vtkIdType> pts(N);
 827   for (int i = 0; i < N; ++i) {
 828     points->SetPoint(i, x[i][0], x[i][1], x[i][2]);
 829     pts[i] = i;
 830   }
 831   poly_data->Allocate(1);
 832   poly_data->SetPoints(points);
 833   poly_data->InsertNextCell(VTK_POLYGON, N, pts.data());
 834 }
 835
 836 template <class C>
 837 double EgVtkObject::convexRatio(const C &x, vec3_t n_plane, bool closed_loop)
 838 {
 839   double L_max = -1e99;
 840   double L_min =  1e99;
 841   int N = x.size();
 842   if (closed_loop) {
 843     --N;
 844   }
 845   for (int i = 0; i < N; ++i) {
 846     for (int j = 0; j < N; ++j) {
 847       int p1 = j;
 848       int p2 = j+1;
 849       if (j == N - 1 && !closed_loop) {
 850         p2 = 0;
 851       }
 852       if (i != p1 && i != p2) {
 853         vec3_t n = n_plane.cross(x[p2] - x[p1]);
 854         n.normalise();
 855         double L = (x[i] - x[j])*n;
 856         L_max = max(L, L_max);
 857         L_min = min(L, L_min);
 858       }
 859     }
 860   }
 861   return L_min/L_max;
 862 }
 863
 864 template <typename C>
 865 void EgVtkObject::sharedNodesOfCells(vtkUnstructuredGrid* grid, vtkIdType id_cell1, vtkIdType id_cell2, C& cont)
 866 {
 867   cont.clear();
 868   EG_GET_CELL(id_cell1, grid);
 869   for (int i = 0; i < num_pts; ++i) {
 870     if (cellContainsNode(grid, id_cell2, pts[i])) {
 871       cont << pts[i];
 872     }
 873   }
 874 }
 875
 876 template <typename C>
 877 void EgVtkObject::getPointsOfCell(vtkUnstructuredGrid* grid, vtkIdType id_cell, C& cont)
 878 {
 879   cont.clear();
 880   EG_VTKSP(vtkIdList, stream);
 881   vtkIdType type_cell = grid->GetCellType(id_cell);
 882   if (type_cell == VTK_POLYHEDRON) {
 883     grid->GetFaceStream(id_cell, stream);
 884     vtkIdType id = 1;
 885     for (int i = 0; i < stream->GetId(0); ++i) {
 886       int num_pts = stream->GetId(id++);
 887       for (int j = 0; j < num_pts; ++j) {
 888         cont << stream->GetId(id);
 889         ++id;
 890       }
 891     }
 892   } else {
 893     grid->GetCellPoints(id_cell, stream);
 894     for (vtkIdType i = 0; i < stream->GetNumberOfIds(); ++i) {
 895       cont << stream->GetId(i);
 896     }
 897   }
 898 }
 899
 900 template <typename C>
 901 void EgVtkObject::vtkIdList2QContainer(vtkIdList *id_list, C &cont)
 902 {
 903   cont.clear();
 904   for (int i = 0; i < id_list->GetNumberOfIds(); ++i) {
 905     cont << id_list->GetId(i);
 906   }
 907 }
 908
 909 template <typename C1, typename C2>
 910 void EgVtkObject::triangulatePolygon(vtkUnstructuredGrid *grid, const C1 &polygon, C2 &triangles)
 911 {
 912   int N = polygon.size();
 913   if (N < 3) {
 914     EG_BUG;
 915   }
 916   QVector<vtkIdType> poly(N+4);
 917   for (int i = 2; i <= N+1; ++i) {
 918     poly[i] = polygon[i-2];
 919   }
 920   poly[0]   = poly[N];
 921   poly[1]   = poly[N+1];
 922   poly[N+2] = poly[2];
 923   poly[N+3] = poly[3];
 924
 925   int    i_best     = 0;
 926   double best_angle = M_PI;
 927
 928   for (int i = 2; i <= N+1; ++i) {
 929     double angle = 0;
 930     vec3_t a, b, c;
 931     grid->GetPoint(poly[i], a.data());
 932     for (int j = 0; j <= N+1; ++j) {
 933       if (j < i-1 || j > i) {
 934         grid->GetPoint(poly[j],   b.data());
 935         grid->GetPoint(poly[j+1], c.data());
 936         angle = max(angle, GeometryTools::angle(b-a, c-a));
 937         angle = max(angle, GeometryTools::angle(a-b, c-b));
 938         angle = max(angle, GeometryTools::angle(a-c, b-c));
 939       }
 940     }
 941     if (angle < best_angle) {
 942       i_best = i;
 943       best_angle = angle;
 944     }
 945   }
 946   triangles.resize(N-2);
 947   int i1 = 0;
 948   int i2 = 1;
 949   if (i1 == i_best || i1 == i_best - N) {
 950     ++i1;
 951     ++i2;
 952   }
 953   for (int i = 0; i < triangles.size(); ++i) {
 954     triangles[i].resize(3);
 955     if (i2 == i_best || i2 == i_best - N) {
 956       i1 += 2;
 957       i2 += 2;
 958     }
 959     triangles[i][0] = poly[i1];
 960     triangles[i][1] = poly[i2];
 961     triangles[i][2] = poly[i_best];
 962     ++i1;
 963     ++i2;
 964   }
 965 }
 966
 967 template <typename C>
 968 void EgVtkObject::invertQContainer(C &cont)
 969 {
 970   QList<typename C::value_type> original;
 971   foreach (typename C::value_type item, cont) {
 972     original << item;
 973   }
 974   cont.clear();
 975   while (original.size() > 0) {
 976     cont << original.last();
 977     original.pop_back();
 978   }
 979 }
 980
 981 template <typename C>
 982 void EgVtkObject::getSurfaceCells(const C &bcs,  QVector<vtkIdType>  &cells,  vtkUnstructuredGrid *grid)
 983 {
 984   int N = 0;
 985   EG_VTKDCC(vtkIntArray, cell_code, grid, "cell_code");
 986   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 987     if (isSurface(id_cell, grid)) {
 988       if (bcs.contains(cell_code->GetValue(id_cell))) {
 989         ++N;
 990       }
 991     }
 992   }
 993   cells.resize(N);
 994   N = 0;
 995   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 996     if (isSurface(id_cell, grid)) {
 997       if (bcs.contains(cell_code->GetValue(id_cell))) {
 998         cells[N] = id_cell;
 999         ++N;
1000       }
1001     }
1002   }
1003 }
1004
1005 template <>
1006 inline bool EgVtkObject::getXmlSetting<QString>(QString key, QString xml_section, QString &value)
1007 {
1008   QString buffer = getXmlSection(xml_section);
1009   bool found = false;
1010   QStringList items = buffer.split(";", QString::SkipEmptyParts);
1011   foreach (QString item, items) {
1012     QStringList words = item.split("=", QString::SkipEmptyParts);
1013     if (words.size() == 2) {
1014       if (words[0] == key) {
1015         found = true;
1016         value = words[1];
1017         break;
1018       }
1019     }
1020   }
1021   return found;
1022 }
1023
1024 template <>
1025 inline bool EgVtkObject::getXmlSetting<double>(QString key, QString xml_section, double &value)
1026 {
1027   QString value_text;
1028   bool found = getXmlSetting(key, xml_section, value_text);
1029   value = value_text.toDouble();
1030   return found;
1031 }
1032
1033 template <>
1034 inline bool EgVtkObject::getXmlSetting<float>(QString key, QString xml_section, float &value)
1035 {
1036   QString value_text;
1037   bool found = getXmlSetting(key, xml_section, value_text);
1038   value = value_text.toFloat();
1039   return found;
1040 }
1041
1042 #endif