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    * Get all surface cells of a grid with a specific boundary condition.
 283    * @param bc    The boundary conditions
 284    * @param cells On return this will hold the Ids of the surface cells.
 285    * @param grid  The grid to operate on.
 286    */
 287   void getSurfaceCells(int bc, QVector<vtkIdType> &cells, vtkUnstructuredGrid *grid);
 288
 289   /**
 290    * Create a cell neighbourship list for a subset grid.
 291    * This has been implemented using VTK's vtkCellLinks structures.
 292    * @param cells the subset of cells
 293    * @param c2c   On return this will hold the neighbourship list
 294    * @param grid  The grid to operate on.
 295    */
 296   void createCellToCell(QVector<vtkIdType> &cells, QVector<QVector<int> > &c2c, vtkUnstructuredGrid *grid);
 297
 298   /**
 299    * Insert a subset of a grid into a vtkPolyData structure.
 300    * This is can be used in order to make use of many readily available
 301    * operations within VTK; one example is smoothing.
 302    * Cell index and node index arrays will be created and passed to the
 303    * poly data structure. Thus any purely geometric effect (no topology change)
 304    * can be directly reintroduced into the vtkUnstructuredGrid.
 305    * @param cells the subset of cells
 306    * @param pdata the vtkPolyData to add the nodes and cells to
 307    * @param grid  The grid to operate on.
 308    */
 309   void addToPolyData(QVector<vtkIdType>  &cells, vtkPolyData *pdata, vtkUnstructuredGrid *grid);
 310
 311   /**
 312    * Copy the attributes from an input to an output cell.
 313    * @param old_grid the input grid
 314    * @param oldId    the existing input cell
 315    * @param new_grid the output grid
 316    * @param newId    the new output cell
 317    */
 318   void copyCellData(vtkUnstructuredGrid *old_grid, vtkIdType oldId, vtkUnstructuredGrid *new_grid, vtkIdType newId);
 319
 320   /**
 321    * Copy the attributes from an input to an output node.
 322    * @param old_grid the input grid
 323    * @param oldId    the existing input node
 324    * @param new_grid the output grid
 325    * @param newId    the new output node
 326    */
 327   void copyNodeData(vtkUnstructuredGrid *old_grid, vtkIdType oldId, vtkUnstructuredGrid *new_grid, vtkIdType newId);
 328
 329   /**
 330    * Create the basic fields on a given grid.
 331    * Care should be taken with the overwrite parameter; if it is set to <i>false</i>
 332    * and the field does not have the correct size it can lead to <i>segmentation faults</i>.
 333    * @param Ncells the number of output cells
 334    * @param Nnodes the number of output nodes
 335    * @param overwrite f set to true existing fields will be re-created
 336    */
 337   void createBasicFields(vtkUnstructuredGrid *grid, vtkIdType num_cells = -1, vtkIdType num_nodes = -1, bool overwrite = false);
 338
 339   /**
 340    * Create the basic cell fields on a given grid.
 341    * Care should be taken with the overwrite parameter; if it is set to <i>false</i>
 342    * and the field does not have the correct size it can lead to <i>segmentation faults</i>.
 343    * @param Ncells the number of output cells
 344    * @param overwrite f set to true existing fields will be re-created
 345    */
 346   void createBasicCellFields(vtkUnstructuredGrid *grid, vtkIdType Ncells, bool overwrite = false);
 347
 348   /**
 349    * Create the basic node fields on a given grid.
 350    * Care should be taken with the overwrite parameter; if it is set to <i>false</i>
 351    * and the field does not have the correct size it can lead to <i>segmentation faults</i>.
 352    * @param Nnodes the number of output nodes
 353    * @param overwrite f set to true existing fields will be re-created
 354    */
 355   void createBasicNodeFields(vtkUnstructuredGrid *grid, vtkIdType Nnodes, bool overwrite = false);
 356
 357   /**
 358    * Allocate memory for a grid. This method will also create the basic
 359    * attribute fields (e.g. "cell_code").
 360    * @param grid   the grid for which to allocate memory
 361    * @param Ncells the number of output cells
 362    * @param Nnodes the number of output nodes
 363    * @param create_fields flag to determine if node and cell data shall be created
 364    */
 365   void allocateGrid(vtkUnstructuredGrid *grid, vtkIdType Ncells, vtkIdType Nnodes, bool create_fields = true);
 366
 367   /**
 368    * Get the names of all node (point) attributes (fields) of a VTK grid
 369    * @param field_names On return this vector will contain the names of all fields
 370    * @param grid the grid
 371    */
 372   void getAllNodeDataNames(QVector<QString> &field_names, vtkUnstructuredGrid *grid);
 373
 374   /**
 375    * Get the names of all cell attributes (fields) of a VTK grid
 376    * @param field_names On return this vector will contain the names of all fields
 377    * @param grid the grid
 378    */
 379   void getAllCellDataNames(QVector<QString> &field_names, vtkUnstructuredGrid *grid);
 380
 381   /**
 382    * Compute the intersection of two Q containers.
 383    * This will return a set.
 384    * @param set1 the first container
 385    * @param set2 the second container
 386    * @param inters on return this will hold the intersection
 387    */
 388   template <class C1, class C2>
 389   void qcontIntersection(const C1& c1, const C2& c2, QSet<typename C1::value_type> &inters);
 390
 391   /**
 392    * Compute the intersection of two Q containers.
 393    * This will return a vector.
 394    * @param set1 the first container
 395    * @param set2 the second container
 396    * @param inters on return this will hold the intersection
 397    */
 398   template <class C1, class C2>
 399   void qcontIntersection(const C1& c1, const C2& c2, QVector<typename C1::value_type> &inters);
 400
 401   /**
 402    * Compute the centre of a cell
 403    * @param grid the grid to use
 404    * @param id_cell the cell of which to compute the centre
 405    * @return the centre of the cell
 406    */
 407   vec3_t cellCentre(vtkUnstructuredGrid *grid, vtkIdType id_cell);
 408
 409   /**
 410    * Compute the angle between two faces.
 411    * @param grid the grid to use
 412    * @param id_face1 index of the first face
 413    * @param id_face2 index of the second face
 414    * @return the angle between the faces (M_PI for a flat surface)
 415    */
 416   double faceAngle(vtkUnstructuredGrid* grid, vtkIdType id_face1, vtkIdType id_face2);
 417
 418   /**
 419    * Get the cells of a grid that are not part of a given set of cells.
 420    * @param grid the grid to use
 421    * @param cells the given set of cells
 422    * @param rest_cells on return this will hold the rest of all cells of the grid (not part of cells)
 423    */
 424   void getRestCells(vtkUnstructuredGrid *grid, const QVector<vtkIdType> &cells, QVector<vtkIdType> &rest_cells);
 425
 426   /**
 427    * Find the corresponding volume cell of a surface cell
 428    * @param grid the grid to use
 429    * @param id_surf the id of the surface cell
 430    * @param n2n the node to cell structure for this grid
 431    * @return the id of the corresponding volume cell (or -1 if not found)
 432    */
 433   vtkIdType findVolumeCell(vtkUnstructuredGrid *grid, vtkIdType id_surf, g2l_t _nodes, l2g_t cells, g2l_t _cells, l2l_t n2c);
 434
 435   /**
 436    * @brief Copy a cell from one grid to another.
 437    * @param src the source grid
 438    * @param id_cell the cell index within the source grid
 439    * @param dst the destination grid
 440    * @return the index of the cell in the destination grid
 441    */
 442   vtkIdType copyCell(vtkUnstructuredGrid* src, vtkIdType id_cell, vtkUnstructuredGrid* dst, vtkIdType offset = 0);
 443
 444   /**
 445    * @brief Copy a cell from one grid to another and translate node indices.
 446    * @param src the source grid
 447    * @param id_cell the cell index within the source grid
 448    * @param dst the destination grid
 449    * @param src2dst a generic container containing the mapping from src to destination index
 450    * @return the index of the cell in the destination grid
 451    */
 452   template <typename C>
 453   vtkIdType copyCell(vtkUnstructuredGrid* src, vtkIdType id_cell, vtkUnstructuredGrid* dst, const C& src2dst);
 454
 455   /**
 456    * Copy "src" grid to "dst" grid. Allocate "dst" so that it fits the data of "src".
 457    * @param src a pointer to the source grid
 458    * @param dst a pointer to the destination grid
 459    */
 460   void makeCopy(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst, bool copy_data = true);
 461
 462   /**
 463    * Copy a part of "src" grid to "dst" grid. Allocate "dst" so that it fits the data to be copied.
 464    * @param src a pointer to the source grid
 465    * @param dst a pointer to the destination grid
 466    * @param cells a container with the cells to be copied
 467    */
 468   template <class C>
 469   void makeCopy(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst, const C &cells);
 470
 471   /**
 472    * Copy "src" grid to "dst" grid. DO NOT allocate "dst" so that it fits the data of "src".
 473    * Allocation is left for the user to do.
 474    * @param src a pointer to the source grid
 475    * @param dst a pointer to the destination grid
 476    */
 477   void makeCopyNoAlloc(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst);
 478
 479   /**
 480    * Change the orientation of a face.
 481    * @param grid the grid to use
 482    * @param id_face the id of the face to change
 483    */
 484   void reorientateFace(vtkUnstructuredGrid *grid, vtkIdType id_face);
 485
 486   /**
 487    * Reset face orientation to original orientation.
 488    * @param grid the grid with the faces
 489    */
 490   void resetOrientation(vtkUnstructuredGrid *grid);
 491
 492   void createIndices(vtkUnstructuredGrid *grid);
 493
 494   /**
 495    * Get the boundary condition of a boundary code.
 496    * @param bc the boundary code
 497    * @return the boundary condition
 498    */
 499   BoundaryCondition getBC(int bc);
 500
 501   /**
 502    * Save the subgrid defined by cls from grid.
 503    * @param grid The source grid
 504    * @param cls The cells to extract
 505    * @param file_name The file to save to
 506    */
 507   template <class C>
 508   void writeCells(vtkUnstructuredGrid *grid, const C &cls, QString file_name);
 509
 510   /**
 511    * Get the SubGrid defined by cls from grid. The function takes care of allocation for SubGrid.
 512    * @param grid The source grid
 513    * @param cls The cells to extract
 514    * @param SubGrid The SubGrid to create
 515    */
 516   template <class C>
 517   void getSubGrid(vtkUnstructuredGrid *grid, const C &cls, vtkUnstructuredGrid *SubGrid);
 518
 519
 520   /**
 521    * Save grid to file filename.
 522    * @param grid The source grid
 523    * @param filename Name of the file to save to
 524    */
 525   void writeGrid(vtkUnstructuredGrid *grid, QString filename);
 526
 527   /**
 528    * Get a file name without extension.
 529    * @param file_name the full name (with extension)
 530    * @return the name without the extension
 531    */
 532   QString stripFromExtension(QString file_name);
 533
 534   /**
 535    * Get the extension of a file name
 536    * @param file_name the full name (with extension)
 537    * @return the extension
 538    */
 539   QString getExtension(QString file_name);
 540
 541   /**
 542    * Get a face of a cell
 543    * @param grid the unstructured grid
 544    * @param id_cell the index of the cell
 545    * @param i_face the index of the face within the cell
 546    * @param ids on return this will contain the nodes of the face
 547    */
 548   void getFaceOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_face, QVector<vtkIdType> &ids);
 549
 550   /**
 551    * Get the normal of a face of a volume cell.
 552    * @param grid the unstructured grid
 553    * @param id_cell the index of the cell
 554    * @param i_face the index of the face within the cell
 555    * @return the normal vector (absolute value corresponds to the area)
 556    */
 557   vec3_t getNormalOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_face);
 558
 559   /**
 560    * Get the centre of a face of a volume cell.
 561    * @param grid the unstructured grid
 562    * @param id_cell the index of the cell
 563    * @param i_face the index of the face within the cell
 564    * @return the centre
 565    */
 566   vec3_t getCentreOfCellFace(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_face);
 567
 568   /**
 569    * Get an edge of a face/cell
 570    * @param grid the unstructured grid
 571    * @param id_cell the index of the cell
 572    * @param i_edge the index of the edge within the cell
 573    * @param ids on return this will contain the nodes of the edge
 574    */
 575   void getEdgeOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_edge, QVector<vtkIdType> &ids);
 576
 577   /**
 578    * Get all boundary codes fo a grid.
 579    * @param grid the grid to extract the boundaru codes from
 580    * @return a set with all boundary codes
 581    */
 582   QSet<int> getAllBoundaryCodes(vtkUnstructuredGrid *grid);
 583
 584   /**
 585    * Check if a cell (face) contains a given node
 586    * @param grid the unstructured grid
 587    * @param id_cell the id of the cell to investigate
 588    * @param id_node the id of the required node
 589    * @return true if id_cell contains id_node
 590    */
 591   bool cellContainsNode(vtkUnstructuredGrid *grid, vtkIdType id_cell, vtkIdType id_node);
 592
 593   /**
 594    * Get all node indices which are shared by two cells.
 595    * These cells can be surface or volume cells; also a combination
 596    * of a volume and a surface cell is possible.
 597    * @param grid the grid to use
 598    * @param id_cell1 index of the first cell
 599    * @param id_cell2 index of the second cell
 600    * @param cont a generic Qt container which will hold the shared node indices on return
 601    */
 602   template <typename C>
 603   void sharedNodesOfCells(vtkUnstructuredGrid* grid, vtkIdType id_cell1, vtkIdType id_cell2, C& cont);
 604
 605   /**
 606    * @brief Get all nodes of a cell from a vtkUnstructuredGrid
 607    * This methods collects all nodes of a cell in a generic Qt container.
 608    * It can be used uniformly for VTK_POLYHEDRON cells and standard cells.
 609    * @param grid the grid to use
 610    * @param id_cell the cell index
 611    * @param cont a generic Qt container which will hold node indices on return
 612    */
 613   template <typename C>
 614   void getPointsOfCell(vtkUnstructuredGrid* grid, vtkIdType id_cell, C& cont);
 615
 616   template <class C> void createPolyData(const C &x, vtkPolyData *poly_data, bool closed_loop = false);
 617   void createPolyDataC2C(vtkPolyData *poly_data, QVector<QVector<vtkIdType> > &c2c);
 618   void createPolyDataN2C(vtkPolyData *poly_data, QVector<QSet<vtkIdType> > &n2c);
 619   void createPolyDataN2N(vtkPolyData *poly_data, QVector<QSet<vtkIdType> > &n2n);
 620   template <class C> double convexRatio(const C &x, vec3_t n_plane, bool closed_loop = false);
 621
 622   template <class C> void invertQContainer(C &cont);
 623
 624
 625 public: // methods
 626
 627   EgVtkObject() { DebugLevel = 0; }
 628
 629   void setBoundaryCodes(const QSet<int> &bcs);
 630   QSet<int> getBoundaryCodes();
 631   void setDebugLevel(int a_DebugLevel) { DebugLevel = a_DebugLevel; }
 632
 633   bool saveGrid( vtkUnstructuredGrid* a_grid, QString file_name );
 634
 635   vtkIdType addGrid(vtkUnstructuredGrid *main_grid, vtkUnstructuredGrid *grid_to_add, vtkIdType offset);
 636
 637   template <typename C>
 638   void vtkIdList2QContainer(vtkIdList *id_list, C &cont);
 639
 640   void checkGridConsitency(vtkUnstructuredGrid *grid);
 641
 642   template <typename C1, typename C2>
 643   void triangulatePolygon(vtkUnstructuredGrid *grid, const C1 &polygon, C2 &triangles);
 644
 645 private:
 646
 647   void addVtkTypeInfo(vtkUnstructuredGrid* a_grid); ///< Add VTK type information to the grid (useful for visualisation with ParaView).
 648
 649 };
 650
 651 //End of class EgVtkObject
 652
 653 template <class C1, class C2>
 654 void EgVtkObject::qcontIntersection(const C1& c1, const C2& c2, QSet<typename C1::value_type> &inters)
 655 {
 656   inters.clear();
 657   foreach (typename C1::value_type t1, c1) {
 658     foreach (typename C2::value_type t2, c2) {
 659       if (t1 == t2) {
 660         inters.insert(t1);
 661       }
 662     }
 663   }
 664 }
 665
 666 template <class C1, class C2>
 667 void EgVtkObject::qcontIntersection(const C1& c1, const C2& c2, QVector<typename C1::value_type> &inters)
 668 {
 669   QSet<typename C1::value_type> inters_set;
 670   qcontIntersection(c1, c2, inters_set);
 671   inters.resize(inters_set.size());
 672   qCopy(inters_set.begin(), inters_set.end(), inters.begin());
 673 }
 674
 675 template <class C>
 676 void EgVtkObject::getSubGrid(vtkUnstructuredGrid *grid, const C &cls, vtkUnstructuredGrid *sub_grid)
 677 {
 678   createIndices(grid);
 679   QVector<vtkIdType> cells;
 680   QVector<vtkIdType> nodes;
 681   cells.resize(cls.size());
 682   qCopy(cls.begin(), cls.end(), cells.begin());
 683   getNodesFromCells(cells, nodes, grid);
 684   allocateGrid(sub_grid, cells.size(), nodes.size());
 685   vtkIdType id_new_node = 0;
 686   QVector<vtkIdType> old2new(grid->GetNumberOfPoints(), -1);
 687   foreach (vtkIdType id_node, nodes) {
 688     vec3_t x;
 689     grid->GetPoint(id_node, x.data());
 690     sub_grid->GetPoints()->SetPoint(id_new_node, x.data());
 691     old2new[id_node] = id_new_node;
 692     copyNodeData(grid, id_node, sub_grid, id_new_node);
 693     ++id_new_node;
 694   }
 695   foreach (vtkIdType id_cell, cells) {
 696     vtkIdType id_new_cell = copyCell(grid, id_cell, sub_grid, old2new);
 697     copyCellData(grid, id_cell, sub_grid, id_new_cell);
 698   }
 699 }
 700
 701 template <class C>
 702 void EgVtkObject::writeCells(vtkUnstructuredGrid *grid, const C &cls, QString file_name)
 703 {
 704   qDebug()<<"Saving cells from grid as "<<file_name;
 705
 706   EG_VTKSP(vtkUnstructuredGrid,SubGrid);
 707   getSubGrid(grid,cls,SubGrid);
 708
 709   EG_VTKSP(vtkXMLUnstructuredGridWriter,vtu);
 710   vtu->SetFileName(qPrintable(file_name));
 711   vtu->SetDataModeToBinary();
 712   vtu->SetInputData(SubGrid);
 713   vtu->Write();
 714 }
 715
 716 template <class C>
 717 void EgVtkObject::getNodesFromCells(const C& cells, QVector<vtkIdType>  &nodes, vtkUnstructuredGrid *grid)
 718 {
 719   QSet<vtkIdType> ex_nodes;
 720   vtkIdType id_cell;
 721   foreach(id_cell, cells) {
 722     QList<vtkIdType> pts;
 723     getPointsOfCell(grid, id_cell, pts);
 724     foreach (vtkIdType id_node, pts) {
 725       if (id_node >= grid->GetNumberOfPoints()) {
 726         EG_BUG;
 727       }
 728       ex_nodes.insert(id_node);
 729     }
 730   }
 731   nodes.resize(ex_nodes.size());
 732   {
 733     int j = 0;
 734     vtkIdType i;
 735     foreach(i,ex_nodes) {
 736       nodes[j] = i;
 737       ++j;
 738     }
 739   }
 740 }
 741
 742 inline vtkIdType EgVtkObject::copyCell(vtkUnstructuredGrid *src, vtkIdType id_cell, vtkUnstructuredGrid *dst, vtkIdType offset)
 743 {
 744   EG_VTKSP(vtkIdList, stream);
 745   vtkIdType type_cell = src->GetCellType(id_cell);
 746   vtkIdType id_new_cell = -1;
 747   if (type_cell == VTK_POLYHEDRON) {
 748     src->GetFaceStream(id_cell, stream);
 749     int id = 1;
 750     for (int i = 0; i < stream->GetId(0); ++i) {
 751       int num_pts = stream->GetId(id);
 752       ++id;
 753       for (int j = 0; j < num_pts; ++j) {
 754         stream->SetId(id, stream->GetId(id) + offset);
 755         ++id;
 756       }
 757     }
 758     id_new_cell = dst->InsertNextCell(type_cell, stream);
 759   } else {
 760     src->GetCellPoints(id_cell, stream);
 761     for (int i = 0; i < stream->GetNumberOfIds(); ++i) {
 762       stream->SetId(i, stream->GetId(i) + offset);
 763     }
 764     id_new_cell = dst->InsertNextCell(type_cell, stream);
 765   }
 766   return id_new_cell;
 767 }
 768
 769 template <typename C>
 770 vtkIdType EgVtkObject::copyCell(vtkUnstructuredGrid *src, vtkIdType id_cell, vtkUnstructuredGrid *dst, const C& src2dst)
 771 {
 772   EG_VTKSP(vtkIdList, stream);
 773   vtkIdType type_cell = src->GetCellType(id_cell);
 774   vtkIdType id_new_cell = -1;
 775   if (type_cell == VTK_POLYHEDRON) {
 776     src->GetFaceStream(id_cell, stream);
 777     int id = 1;
 778     for (int i = 0; i < stream->GetId(0); ++i) {
 779       int num_pts = stream->GetId(id);
 780       ++id;
 781       for (int j = 0; j < num_pts; ++j) {
 782         stream->SetId(id, src2dst[stream->GetId(id)]);
 783         ++id;
 784       }
 785     }
 786     id_new_cell = dst->InsertNextCell(type_cell, stream);
 787   } else {
 788     src->GetCellPoints(id_cell, stream);
 789     for (int i = 0; i < stream->GetNumberOfIds(); ++i) {
 790       if (src2dst[stream->GetId(i)] < 0) {
 791         EG_BUG;
 792       }
 793       if (src2dst[stream->GetId(i)] >= dst->GetNumberOfPoints()) {
 794         EG_BUG;
 795       }
 796       stream->SetId(i, src2dst[stream->GetId(i)]);
 797     }
 798     id_new_cell = dst->InsertNextCell(type_cell, stream);
 799   }
 800   return id_new_cell;
 801 }
 802
 803 template <class C>
 804 void EgVtkObject::makeCopy(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst, const C& cells)
 805 {
 806   QVector<vtkIdType> nodes;
 807   getNodesFromCells(cells, nodes, src);
 808   allocateGrid(dst, cells.size(), nodes.size());
 809   vtkIdType id_new_node = 0;
 810   QVector<vtkIdType> old2new(src->GetNumberOfPoints(), -1);
 811   foreach (vtkIdType id_node, nodes) {
 812     vec3_t x;
 813     src->GetPoints()->GetPoint(id_node, x.data());
 814     dst->GetPoints()->SetPoint(id_new_node, x.data());
 815     copyNodeData(src, id_node, dst, id_new_node);
 816     old2new[id_node] = id_new_node;
 817     ++id_new_node;
 818   }
 819   foreach (vtkIdType id_cell, cells) {
 820     vtkIdType id_new_cell = copyCell(src, id_cell, dst, old2new);
 821     copyCellData(src, id_cell, dst, id_new_cell);
 822   }
 823 }
 824
 825 template <class C>
 826 void EgVtkObject::createPolyData(const C &x, vtkPolyData *poly_data, bool closed_loop)
 827 {
 828   int N = x.size();
 829   if (closed_loop) {
 830     --N;
 831   }
 832   EG_VTKSP(vtkPoints, points);
 833   points->SetNumberOfPoints(N);
 834   QVector<vtkIdType> pts(N);
 835   for (int i = 0; i < N; ++i) {
 836     points->SetPoint(i, x[i][0], x[i][1], x[i][2]);
 837     pts[i] = i;
 838   }
 839   poly_data->Allocate(1);
 840   poly_data->SetPoints(points);
 841   poly_data->InsertNextCell(VTK_POLYGON, N, pts.data());
 842 }
 843
 844 template <class C>
 845 double EgVtkObject::convexRatio(const C &x, vec3_t n_plane, bool closed_loop)
 846 {
 847   double L_max = -1e99;
 848   double L_min =  1e99;
 849   int N = x.size();
 850   if (closed_loop) {
 851     --N;
 852   }
 853   for (int i = 0; i < N; ++i) {
 854     for (int j = 0; j < N; ++j) {
 855       int p1 = j;
 856       int p2 = j+1;
 857       if (j == N - 1 && !closed_loop) {
 858         p2 = 0;
 859       }
 860       if (i != p1 && i != p2) {
 861         vec3_t n = n_plane.cross(x[p2] - x[p1]);
 862         n.normalise();
 863         double L = (x[i] - x[j])*n;
 864         L_max = max(L, L_max);
 865         L_min = min(L, L_min);
 866       }
 867     }
 868   }
 869   return L_min/L_max;
 870 }
 871
 872 template <typename C>
 873 void EgVtkObject::sharedNodesOfCells(vtkUnstructuredGrid* grid, vtkIdType id_cell1, vtkIdType id_cell2, C& cont)
 874 {
 875   cont.clear();
 876   EG_GET_CELL(id_cell1, grid);
 877   for (int i = 0; i < num_pts; ++i) {
 878     if (cellContainsNode(grid, id_cell2, pts[i])) {
 879       cont << pts[i];
 880     }
 881   }
 882 }
 883
 884 template <typename C>
 885 void EgVtkObject::getPointsOfCell(vtkUnstructuredGrid* grid, vtkIdType id_cell, C& cont)
 886 {
 887   cont.clear();
 888   EG_VTKSP(vtkIdList, stream);
 889   vtkIdType type_cell = grid->GetCellType(id_cell);
 890   if (type_cell == VTK_POLYHEDRON) {
 891     grid->GetFaceStream(id_cell, stream);
 892     vtkIdType id = 1;
 893     for (int i = 0; i < stream->GetId(0); ++i) {
 894       int num_pts = stream->GetId(id++);
 895       for (int j = 0; j < num_pts; ++j) {
 896         cont << stream->GetId(id);
 897         ++id;
 898       }
 899     }
 900   } else {
 901     grid->GetCellPoints(id_cell, stream);
 902     for (vtkIdType i = 0; i < stream->GetNumberOfIds(); ++i) {
 903       cont << stream->GetId(i);
 904     }
 905   }
 906 }
 907
 908 template <typename C>
 909 void EgVtkObject::vtkIdList2QContainer(vtkIdList *id_list, C &cont)
 910 {
 911   cont.clear();
 912   for (int i = 0; i < id_list->GetNumberOfIds(); ++i) {
 913     cont << id_list->GetId(i);
 914   }
 915 }
 916
 917 template <typename C1, typename C2>
 918 void EgVtkObject::triangulatePolygon(vtkUnstructuredGrid *grid, const C1 &polygon, C2 &triangles)
 919 {
 920   int N = polygon.size();
 921   if (N < 3) {
 922     EG_BUG;
 923   }
 924   QVector<vtkIdType> poly(N+4);
 925   for (int i = 2; i <= N+1; ++i) {
 926     poly[i] = polygon[i-2];
 927   }
 928   poly[0]   = poly[N];
 929   poly[1]   = poly[N+1];
 930   poly[N+2] = poly[2];
 931   poly[N+3] = poly[3];
 932
 933   int    i_best     = 0;
 934   double best_angle = M_PI;
 935
 936   for (int i = 2; i <= N+1; ++i) {
 937     double angle = 0;
 938     vec3_t a, b, c;
 939     grid->GetPoint(poly[i], a.data());
 940     for (int j = 0; j <= N+1; ++j) {
 941       if (j < i-1 || j > i) {
 942         grid->GetPoint(poly[j],   b.data());
 943         grid->GetPoint(poly[j+1], c.data());
 944         angle = max(angle, GeometryTools::angle(b-a, c-a));
 945         angle = max(angle, GeometryTools::angle(a-b, c-b));
 946         angle = max(angle, GeometryTools::angle(a-c, b-c));
 947       }
 948     }
 949     if (angle < best_angle) {
 950       i_best = i;
 951       best_angle = angle;
 952     }
 953   }
 954   triangles.resize(N-2);
 955   int i1 = 0;
 956   int i2 = 1;
 957   if (i1 == i_best || i1 == i_best - N) {
 958     ++i1;
 959     ++i2;
 960   }
 961   for (int i = 0; i < triangles.size(); ++i) {
 962     triangles[i].resize(3);
 963     if (i2 == i_best || i2 == i_best - N) {
 964       i1 += 2;
 965       i2 += 2;
 966     }
 967     triangles[i][0] = poly[i1];
 968     triangles[i][1] = poly[i2];
 969     triangles[i][2] = poly[i_best];
 970     ++i1;
 971     ++i2;
 972   }
 973 }
 974
 975 template <typename C>
 976 void EgVtkObject::invertQContainer(C &cont)
 977 {
 978   QList<typename C::value_type> original;
 979   foreach (typename C::value_type item, cont) {
 980     original << item;
 981   }
 982   cont.clear();
 983   while (original.size() > 0) {
 984     cont << original.last();
 985     original.pop_back();
 986   }
 987 }
 988
 989 template <typename C>
 990 void EgVtkObject::getSurfaceCells(const C &bcs,  QVector<vtkIdType>  &cells,  vtkUnstructuredGrid *grid)
 991 {
 992   int N = 0;
 993   EG_VTKDCC(vtkIntArray, cell_code, grid, "cell_code");
 994   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 995     if (isSurface(id_cell, grid)) {
 996       if (bcs.contains(cell_code->GetValue(id_cell))) {
 997         ++N;
 998       }
 999     }
1000   }
1001   cells.resize(N);
1002   N = 0;
1003   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
1004     if (isSurface(id_cell, grid)) {
1005       if (bcs.contains(cell_code->GetValue(id_cell))) {
1006         cells[N] = id_cell;
1007         ++N;
1008       }
1009     }
1010   }
1011 }
1012
1013 template <>
1014 inline bool EgVtkObject::getXmlSetting<QString>(QString key, QString xml_section, QString &value)
1015 {
1016   QString buffer = getXmlSection(xml_section);
1017   bool found = false;
1018   QStringList items = buffer.split(";", QString::SkipEmptyParts);
1019   foreach (QString item, items) {
1020     QStringList words = item.split("=", QString::SkipEmptyParts);
1021     if (words.size() == 2) {
1022       if (words[0] == key) {
1023         found = true;
1024         value = words[1];
1025         break;
1026       }
1027     }
1028   }
1029   return found;
1030 }
1031
1032 template <>
1033 inline bool EgVtkObject::getXmlSetting<double>(QString key, QString xml_section, double &value)
1034 {
1035   QString value_text;
1036   bool found = getXmlSetting(key, xml_section, value_text);
1037   value = value_text.toDouble();
1038   return found;
1039 }
1040
1041 template <>
1042 inline bool EgVtkObject::getXmlSetting<float>(QString key, QString xml_section, float &value)
1043 {
1044   QString value_text;
1045   bool found = getXmlSetting(key, xml_section, value_text);
1046   value = value_text.toFloat();
1047   return found;
1048 }
1049
1050 #endif