src/egvtkobject.cpp

   1 //
   2 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   3 // +                                                                      +
   4 // + This file is part of enGrid.                                         +
   5 // +                                                                      +
   6 // + Copyright 2008,2009 Oliver Gloth                                     +
   7 // +                                                                      +
   8 // + enGrid is free software: you can redistribute it and/or modify       +
   9 // + it under the terms of the GNU General Public License as published by +
  10 // + the Free Software Foundation, either version 3 of the License, or    +
  11 // + (at your option) any later version.                                  +
  12 // +                                                                      +
  13 // + enGrid is distributed in the hope that it will be useful,            +
  14 // + but WITHOUT ANY WARRANTY; without even the implied warranty of       +
  15 // + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the        +
  16 // + GNU General Public License for more details.                         +
  17 // +                                                                      +
  18 // + You should have received a copy of the GNU General Public License    +
  19 // + along with enGrid. If not, see <http://www.gnu.org/licenses/>.       +
  20 // +                                                                      +
  21 // ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  22 //
  23 #include "egvtkobject.h"
  24 #include "guimainwindow.h"
  25 #include "beziertriangle.h"
  26
  27 #include <vtkCellData.h>
  28 #include <vtkPointData.h>
  29 #include <vtkCellLinks.h>
  30 #include <vtkCellType.h>
  31 #include <vtkIdList.h>
  32 #include <vtkCell.h>
  33 #include <vtkCharArray.h>
  34
  35 int EgVtkObject::DebugLevel;
  36
  37 void EgVtkObject::computeNormals
  38 (
  39   QVector<vec3_t>     &cell_normals,
  40   QVector<vec3_t>     &node_normals,
  41   QVector<vtkIdType>  &cells,
  42   QVector<vtkIdType>  &nodes,
  43   vtkUnstructuredGrid *grid
  44 )
  45 {
  46   using namespace GeometryTools;
  47
  48   cell_normals.resize(cells.size());
  49   node_normals.fill(vec3_t(0,0,0), nodes.size());
  50   QVector<int> g2s;
  51   createNodeMapping(nodes, g2s, grid);
  52   for (int i_cell = 0; i_cell < cells.count(); ++i_cell) {
  53     vtkIdType id_cell = cells[i_cell];
  54     vtkIdType *pts;
  55     vtkIdType npts;
  56     grid->GetCellPoints(id_cell, npts, pts);
  57     cell_normals[i_cell] = cellNormal(grid, id_cell);
  58     cell_normals[i_cell].normalise();
  59     for (int i_pts = 0; i_pts < npts; ++i_pts) {
  60       if (g2s[pts[i_pts]] != -1) {
  61         node_normals[g2s[pts[i_pts]]] += cell_normals[i_cell];
  62       }
  63     }
  64   }
  65   for (int i_node = 0; i_node < nodes.count(); ++i_node) {
  66     node_normals[i_node].normalise();
  67     //cout << node_normals[i_node] << endl;
  68   }
  69 }
  70
  71 void EgVtkObject::createNodeMapping
  72 (
  73   QVector<vtkIdType>  &nodes,
  74   QVector<int>        &_nodes,
  75   vtkUnstructuredGrid *grid
  76 )
  77 {
  78   _nodes.fill(-1,grid->GetNumberOfPoints());
  79   for (int i_nodes = 0; i_nodes < nodes.size(); ++i_nodes) {
  80     _nodes[nodes[i_nodes]] = i_nodes;
  81   }
  82 }
  83
  84 void EgVtkObject::createCellMapping
  85 (
  86   QVector<vtkIdType>  &cells,
  87   QVector<int>        &_cells,
  88   vtkUnstructuredGrid *grid
  89 )
  90 {
  91   _cells.fill(-1,grid->GetNumberOfCells());
  92   for (int i_cells = 0; i_cells < cells.size(); ++i_cells) {
  93     _cells[cells[i_cells]] = i_cells;
  94   }
  95 }
  96
  97 void EgVtkObject::createNodeToBcMapping
  98 (
  99   QVector<QSet<int> >  &bcs,
 100   vtkUnstructuredGrid  *grid
 101 )
 102 {
 103   EG_BUG;
 104   bcs.fill(QSet<int>(), grid->GetNumberOfPoints());
 105   grid->BuildLinks();
 106   EG_VTKDCC(vtkIntArray, cell_code, grid, "cell_code");
 107   for (vtkIdType nodeId = 0; nodeId < grid->GetNumberOfPoints(); ++nodeId) {
 108     int Ncells = grid->GetCellLinks()->GetNcells(nodeId);
 109     for (int i = 0; i < Ncells; ++i) {
 110       vtkIdType id_cell = grid->GetCellLinks()->GetCells(nodeId)[i];
 111       vtkIdType ct = grid->GetCellType(id_cell);
 112       if ((ct == VTK_TRIANGLE) || (ct = VTK_QUAD)) {
 113         if (cell_code->GetValue(id_cell) > 0) {
 114           bcs[nodeId].insert(cell_code->GetValue(id_cell));
 115         }
 116       }
 117     }
 118   }
 119 }
 120
 121 void EgVtkObject::createNodeToCell
 122 (
 123   QVector<vtkIdType>  &cells,
 124   QVector<vtkIdType>  &nodes,
 125   QVector<int>        &_nodes,
 126   QVector<QSet<int> > &n2c,
 127   vtkUnstructuredGrid *grid
 128 )
 129 {
 130   n2c.fill(QSet<int>(), nodes.size());
 131   for (vtkIdType i_cells = 0; i_cells < cells.size(); ++i_cells) {
 132     vtkIdType *pts;
 133     vtkIdType  Npts;
 134     grid->GetCellPoints(cells[i_cells], Npts, pts);
 135     for (int i_pts = 0; i_pts < Npts; ++i_pts) {
 136       n2c[_nodes[pts[i_pts]]].insert(i_cells);
 137     }
 138   }
 139 }
 140
 141 void EgVtkObject::createNodeToCell
 142 (
 143   QVector<vtkIdType>     &cells,
 144   QVector<vtkIdType>     &nodes,
 145   QVector<int>           &_nodes,
 146   QVector<QVector<int> > &n2c,
 147   vtkUnstructuredGrid    *grid
 148 )
 149 {
 150   n2c.fill(QVector<int>(), nodes.size());
 151   QVector<int> count(nodes.size(),0);
 152   for (vtkIdType i_cells = 0; i_cells < cells.size(); ++i_cells) {
 153     vtkIdType *pts;
 154     vtkIdType  Npts;
 155     grid->GetCellPoints(cells[i_cells], Npts, pts);
 156     for (int i_pts = 0; i_pts < Npts; ++i_pts) {
 157       ++count[_nodes[pts[i_pts]]];
 158     }
 159   }
 160   for (int i = 0; i < nodes.size(); ++i) {
 161     n2c[i].resize(count[i]);
 162     count[i] = 0;
 163   }
 164   for (vtkIdType i_cells = 0; i_cells < cells.size(); ++i_cells) {
 165     vtkIdType *pts;
 166     vtkIdType  Npts;
 167     grid->GetCellPoints(cells[i_cells], Npts, pts);
 168     for (int i_pts = 0; i_pts < Npts; ++i_pts) {
 169       int i_nodes = _nodes[pts[i_pts]];
 170       n2c[i_nodes][count[i_nodes]] = i_cells;
 171       ++count[i_nodes];
 172     }
 173   }
 174 }
 175
 176 void EgVtkObject::addToN2N(QVector<QSet<int> > &n2n, int n1, int n2)
 177 {
 178   n2n[n1].insert(n2);
 179   n2n[n2].insert(n1);
 180 }
 181
 182 void EgVtkObject::createNodeToNode(QVector<vtkIdType> &cells, QVector<vtkIdType> &nodes, QVector<int> &_nodes, QVector<QSet<int> > &n2n, vtkUnstructuredGrid *grid)
 183 {
 184   n2n.fill(QSet<int>(), nodes.size());
 185   foreach (vtkIdType id_cell, cells) {
 186     vtkIdType *pts;
 187     vtkIdType  Npts;
 188     grid->GetCellPoints(id_cell, Npts, pts);
 189     vector<int> n(Npts);
 190     for (int i = 0; i < Npts; ++i) {
 191       n[i] = _nodes[pts[i]];
 192     }
 193     vtkIdType cellType = grid->GetCellType(id_cell);
 194     if (cellType == VTK_TRIANGLE) {
 195       addToN2N(n2n, n[0], n[1]);
 196       addToN2N(n2n, n[1], n[2]);
 197       addToN2N(n2n, n[2], n[0]);
 198     } else if (cellType == VTK_QUAD) {
 199       addToN2N(n2n, n[0], n[1]);
 200       addToN2N(n2n, n[1], n[2]);
 201       addToN2N(n2n, n[2], n[3]);
 202       addToN2N(n2n, n[3], n[0]);
 203     } else if (cellType == VTK_TETRA) {
 204       addToN2N(n2n, n[0], n[1]);
 205       addToN2N(n2n, n[0], n[2]);
 206       addToN2N(n2n, n[0], n[3]);
 207       addToN2N(n2n, n[1], n[2]);
 208       addToN2N(n2n, n[1], n[3]);
 209       addToN2N(n2n, n[2], n[3]);
 210     } else if (cellType == VTK_PYRAMID) {
 211       addToN2N(n2n, n[0], n[1]);
 212       addToN2N(n2n, n[0], n[3]);
 213       addToN2N(n2n, n[0], n[4]);
 214       addToN2N(n2n, n[1], n[2]);
 215       addToN2N(n2n, n[1], n[4]);
 216       addToN2N(n2n, n[2], n[3]);
 217       addToN2N(n2n, n[2], n[4]);
 218       addToN2N(n2n, n[3], n[4]);
 219     } else if (cellType == VTK_WEDGE) {
 220       addToN2N(n2n, n[0], n[1]);
 221       addToN2N(n2n, n[0], n[2]);
 222       addToN2N(n2n, n[0], n[3]);
 223       addToN2N(n2n, n[1], n[2]);
 224       addToN2N(n2n, n[1], n[4]);
 225       addToN2N(n2n, n[2], n[5]);
 226       addToN2N(n2n, n[3], n[4]);
 227       addToN2N(n2n, n[3], n[5]);
 228       addToN2N(n2n, n[4], n[5]);
 229     } else if (cellType == VTK_HEXAHEDRON) {
 230       addToN2N(n2n, n[0], n[1]);
 231       addToN2N(n2n, n[0], n[3]);
 232       addToN2N(n2n, n[0], n[4]);
 233       addToN2N(n2n, n[1], n[2]);
 234       addToN2N(n2n, n[1], n[5]);
 235       addToN2N(n2n, n[2], n[3]);
 236       addToN2N(n2n, n[2], n[6]);
 237       addToN2N(n2n, n[3], n[7]);
 238       addToN2N(n2n, n[4], n[5]);
 239       addToN2N(n2n, n[4], n[7]);
 240       addToN2N(n2n, n[5], n[6]);
 241       addToN2N(n2n, n[6], n[7]);
 242     }
 243   }
 244 }
 245
 246 void EgVtkObject::createNodeToNode
 247 (
 248   QVector<vtkIdType>     &cells,
 249   QVector<vtkIdType>     &nodes,
 250   QVector<int>           &_nodes,
 251   QVector<QVector<int> > &n2n,
 252   vtkUnstructuredGrid    *grid
 253 )
 254 {
 255   QVector<QSet<int> > n2n_set;
 256   createNodeToNode(cells, nodes, _nodes, n2n_set, grid);
 257   n2n.resize(n2n_set.size());
 258   for (int i = 0; i < n2n.size(); ++i) {
 259     n2n[i].resize(n2n_set[i].size());
 260     qCopy(n2n_set[i].begin(), n2n_set[i].end(), n2n[i].begin());
 261   }
 262 }
 263
 264 void EgVtkObject::getAllCells
 265 (
 266   QVector<vtkIdType>  &cells,
 267   vtkUnstructuredGrid *grid
 268 )
 269 {
 270   int N = 0;
 271   cells.resize(grid->GetNumberOfCells());
 272   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 273     cells[N] = id_cell;
 274     ++N;
 275   }
 276 }
 277
 278 void EgVtkObject::getAllCellsOfType
 279 (
 280   vtkIdType            type,
 281   QVector<vtkIdType>  &cells,
 282   vtkUnstructuredGrid *grid
 283 )
 284 {
 285   int N = 0;
 286   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 287     if (grid->GetCellType(id_cell) == type) {
 288       ++N;
 289     }
 290   }
 291   cells.resize(N);
 292   N = 0;
 293   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 294     if (grid->GetCellType(id_cell) == type) {
 295       cells[N] = id_cell;
 296       ++N;
 297     }
 298   }
 299 }
 300
 301
 302 void EgVtkObject::getAllVolumeCells
 303 (
 304   QVector<vtkIdType>  &cells,
 305   vtkUnstructuredGrid *grid
 306 )
 307 {
 308   int N = 0;
 309   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 310     if (isVolume(id_cell, grid)) {
 311       ++N;
 312     }
 313   }
 314   cells.resize(N);
 315   N = 0;
 316   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 317     if (isVolume(id_cell, grid)) {
 318       cells[N] = id_cell;
 319       ++N;
 320     }
 321   }
 322 }
 323
 324 void EgVtkObject::getAllSurfaceCells
 325 (
 326   QVector<vtkIdType>  &cells,
 327   vtkUnstructuredGrid *grid
 328 )
 329 {
 330   int N = 0;
 331   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 332     if (isSurface(id_cell, grid)) {
 333       ++N;
 334     }
 335   }
 336   cells.resize(N);
 337   N = 0;
 338   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 339     if (isSurface(id_cell, grid)) {
 340       cells[N] = id_cell;
 341       ++N;
 342     }
 343   }
 344 }
 345
 346 void EgVtkObject::getSurfaceCells
 347 (
 348   QSet<int>           &bcs,
 349   QVector<vtkIdType>  &cells,
 350   vtkUnstructuredGrid *grid
 351 )
 352 {
 353   int N = 0;
 354   EG_VTKDCC(vtkIntArray, cell_code, grid, "cell_code");
 355   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 356     if (isSurface(id_cell, grid)) {
 357       if (bcs.contains(cell_code->GetValue(id_cell))) {
 358         ++N;
 359       }
 360     }
 361   }
 362   cells.resize(N);
 363   N = 0;
 364   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 365     if (isSurface(id_cell, grid)) {
 366       if (bcs.contains(cell_code->GetValue(id_cell))) {
 367         cells[N] = id_cell;
 368         ++N;
 369       }
 370     }
 371   }
 372 }
 373
 374 void EgVtkObject::addToC2C(vtkIdType id_cell, QVector<int> &_cells, QVector<QVector<int> > &c2c, int j, vtkIdList *nds, vtkIdList *cls, vtkUnstructuredGrid *grid)
 375 {
 376   c2c[_cells[id_cell]][j] = -1;
 377   grid->GetCellNeighbors(id_cell, nds, cls);
 378   for (int i = 0; i < cls->GetNumberOfIds(); ++i) {
 379     if (cls->GetId(i) != id_cell) {
 380       if (_cells[cls->GetId(i)] != -1) {
 381         c2c[_cells[id_cell]][j] = _cells[cls->GetId(i)];
 382       }
 383     }
 384   }
 385 }
 386
 387
 388 void EgVtkObject::createCellToCell(QVector<vtkIdType> &cells, QVector<QVector<int> > &c2c, vtkUnstructuredGrid *grid)
 389 {
 390   // GetCellNeighbors(vtkIdType id_cell, vtkIdList *ptIds, vtkIdList *id_cells)
 391   grid->BuildLinks();
 392   QVector<int> _cells;
 393   createCellMapping(cells, _cells, grid);
 394   c2c.fill(QVector<int>(), cells.size());
 395   EG_VTKSP(vtkIdList, nds);
 396   EG_VTKSP(vtkIdList, cls);
 397   for (int i = 0; i < cells.size(); ++i) {
 398     vtkIdType id_cell = cells[i];
 399     vtkIdType *pts;
 400     vtkIdType  Npts;
 401     grid->GetCellPoints(id_cell, Npts, pts);
 402     if (grid->GetCellType(id_cell) == VTK_TRIANGLE) {
 403       c2c[i].resize(3);
 404       nds->Reset();
 405       nds->InsertNextId(pts[0]);
 406       nds->InsertNextId(pts[1]);
 407       addToC2C(id_cell, _cells, c2c, 0, nds, cls, grid);
 408       nds->Reset();
 409       nds->InsertNextId(pts[1]);
 410       nds->InsertNextId(pts[2]);
 411       addToC2C(id_cell, _cells, c2c, 1, nds, cls, grid);
 412       nds->Reset();
 413       nds->InsertNextId(pts[2]);
 414       nds->InsertNextId(pts[0]);
 415       addToC2C(id_cell, _cells, c2c, 2, nds, cls, grid);
 416     } else if (grid->GetCellType(id_cell) == VTK_QUAD) {
 417       c2c[i].resize(4);
 418       nds->Reset();
 419       nds->InsertNextId(pts[0]);
 420       nds->InsertNextId(pts[1]);
 421       addToC2C(id_cell, _cells, c2c, 0, nds, cls, grid);
 422       nds->Reset();
 423       nds->InsertNextId(pts[1]);
 424       nds->InsertNextId(pts[2]);
 425       addToC2C(id_cell, _cells, c2c, 1, nds, cls, grid);
 426       nds->Reset();
 427       nds->InsertNextId(pts[2]);
 428       nds->InsertNextId(pts[3]);
 429       addToC2C(id_cell, _cells, c2c, 2, nds, cls, grid);
 430       nds->Reset();
 431       nds->InsertNextId(pts[3]);
 432       nds->InsertNextId(pts[0]);
 433       addToC2C(id_cell, _cells, c2c, 3, nds, cls, grid);
 434     } else if (grid->GetCellType(id_cell) == VTK_TETRA) {
 435       c2c[i].resize(4);
 436       nds->Reset();
 437       nds->InsertNextId(pts[0]);
 438       nds->InsertNextId(pts[1]);
 439       nds->InsertNextId(pts[2]);
 440       addToC2C(id_cell, _cells, c2c, 0, nds, cls, grid);
 441       nds->Reset();
 442       nds->InsertNextId(pts[0]);
 443       nds->InsertNextId(pts[1]);
 444       nds->InsertNextId(pts[3]);
 445       addToC2C(id_cell, _cells, c2c, 1, nds, cls, grid);
 446       nds->Reset();
 447       nds->InsertNextId(pts[0]);
 448       nds->InsertNextId(pts[3]);
 449       nds->InsertNextId(pts[2]);
 450       addToC2C(id_cell, _cells, c2c, 2, nds, cls, grid);
 451       nds->Reset();
 452       nds->InsertNextId(pts[1]);
 453       nds->InsertNextId(pts[2]);
 454       nds->InsertNextId(pts[3]);
 455       addToC2C(id_cell, _cells, c2c, 3, nds, cls, grid);
 456     } else if (grid->GetCellType(id_cell) == VTK_PYRAMID) {
 457       c2c[i].resize(5);
 458       nds->Reset();
 459       nds->InsertNextId(pts[0]);
 460       nds->InsertNextId(pts[1]);
 461       nds->InsertNextId(pts[2]);
 462       nds->InsertNextId(pts[3]);
 463       addToC2C(id_cell, _cells, c2c, 0, nds, cls, grid);
 464       nds->Reset();
 465       nds->InsertNextId(pts[0]);
 466       nds->InsertNextId(pts[1]);
 467       nds->InsertNextId(pts[4]);
 468       addToC2C(id_cell, _cells, c2c, 1, nds, cls, grid);
 469       nds->Reset();
 470       nds->InsertNextId(pts[1]);
 471       nds->InsertNextId(pts[2]);
 472       nds->InsertNextId(pts[4]);
 473       addToC2C(id_cell, _cells, c2c, 2, nds, cls, grid);
 474       nds->Reset();
 475       nds->InsertNextId(pts[2]);
 476       nds->InsertNextId(pts[3]);
 477       nds->InsertNextId(pts[4]);
 478       addToC2C(id_cell, _cells, c2c, 3, nds, cls, grid);
 479       nds->Reset();
 480       nds->InsertNextId(pts[3]);
 481       nds->InsertNextId(pts[0]);
 482       nds->InsertNextId(pts[4]);
 483       addToC2C(id_cell, _cells, c2c, 4, nds, cls, grid);
 484     } else if (grid->GetCellType(id_cell) == VTK_WEDGE) {
 485       c2c[i].resize(5);
 486       nds->Reset();
 487       nds->InsertNextId(pts[0]);
 488       nds->InsertNextId(pts[1]);
 489       nds->InsertNextId(pts[2]);
 490       addToC2C(id_cell, _cells, c2c, 0, nds, cls, grid);
 491       nds->Reset();
 492       nds->InsertNextId(pts[3]);
 493       nds->InsertNextId(pts[4]);
 494       nds->InsertNextId(pts[5]);
 495       addToC2C(id_cell, _cells, c2c, 1, nds, cls, grid);
 496       nds->Reset();
 497       nds->InsertNextId(pts[0]);
 498       nds->InsertNextId(pts[1]);
 499       nds->InsertNextId(pts[4]);
 500       nds->InsertNextId(pts[3]);
 501       addToC2C(id_cell, _cells, c2c, 2, nds, cls, grid);
 502       nds->Reset();
 503       nds->InsertNextId(pts[1]);
 504       nds->InsertNextId(pts[4]);
 505       nds->InsertNextId(pts[5]);
 506       nds->InsertNextId(pts[2]);
 507       addToC2C(id_cell, _cells, c2c, 3, nds, cls, grid);
 508       nds->Reset();
 509       nds->InsertNextId(pts[0]);
 510       nds->InsertNextId(pts[2]);
 511       nds->InsertNextId(pts[5]);
 512       nds->InsertNextId(pts[3]);
 513       addToC2C(id_cell, _cells, c2c, 4, nds, cls, grid);
 514     } else if (grid->GetCellType(id_cell) == VTK_HEXAHEDRON) {
 515       c2c[i].resize(6);
 516       nds->Reset();
 517       nds->InsertNextId(pts[0]);
 518       nds->InsertNextId(pts[3]);
 519       nds->InsertNextId(pts[2]);
 520       nds->InsertNextId(pts[1]);
 521       addToC2C(id_cell, _cells, c2c, 0, nds, cls, grid);
 522       nds->Reset();
 523       nds->InsertNextId(pts[4]);
 524       nds->InsertNextId(pts[5]);
 525       nds->InsertNextId(pts[6]);
 526       nds->InsertNextId(pts[7]);
 527       addToC2C(id_cell, _cells, c2c, 1, nds, cls, grid);
 528       nds->Reset();
 529       nds->InsertNextId(pts[0]);
 530       nds->InsertNextId(pts[1]);
 531       nds->InsertNextId(pts[5]);
 532       nds->InsertNextId(pts[4]);
 533       addToC2C(id_cell, _cells, c2c, 2, nds, cls, grid);
 534       nds->Reset();
 535       nds->InsertNextId(pts[3]);
 536       nds->InsertNextId(pts[7]);
 537       nds->InsertNextId(pts[6]);
 538       nds->InsertNextId(pts[2]);
 539       addToC2C(id_cell, _cells, c2c, 3, nds, cls, grid);
 540       nds->Reset();
 541       nds->InsertNextId(pts[0]);
 542       nds->InsertNextId(pts[4]);
 543       nds->InsertNextId(pts[7]);
 544       nds->InsertNextId(pts[3]);
 545       addToC2C(id_cell, _cells, c2c, 4, nds, cls, grid);
 546       nds->Reset();
 547       nds->InsertNextId(pts[1]);
 548       nds->InsertNextId(pts[2]);
 549       nds->InsertNextId(pts[6]);
 550       nds->InsertNextId(pts[5]);
 551       addToC2C(id_cell, _cells, c2c, 5, nds, cls, grid);
 552     }
 553   }
 554 }
 555
 556 bool EgVtkObject::isVolume(vtkIdType id_cell, vtkUnstructuredGrid *grid)
 557 {
 558   bool isVol = false;
 559   if      (grid->GetCellType(id_cell) == VTK_TETRA)      isVol = true;
 560   else if (grid->GetCellType(id_cell) == VTK_PYRAMID)    isVol = true;
 561   else if (grid->GetCellType(id_cell) == VTK_WEDGE)      isVol = true;
 562   else if (grid->GetCellType(id_cell) == VTK_HEXAHEDRON) isVol = true;
 563   return isVol;
 564 }
 565
 566 bool EgVtkObject::isSurface(vtkIdType id_cell, vtkUnstructuredGrid *grid)
 567 {
 568   bool isSurf = false;
 569   if      (grid->GetCellType(id_cell) == VTK_TRIANGLE) isSurf = true;
 570   else if (grid->GetCellType(id_cell) == VTK_QUAD)     isSurf = true;
 571   return isSurf;
 572 }
 573
 574 void EgVtkObject::UpdateCellIndex(vtkUnstructuredGrid *grid)
 575 {
 576   if (!grid->GetCellData()->GetArray("cell_index")) {
 577     EG_VTKSP(vtkLongArray_t, cell_index);
 578     cell_index->SetName("cell_index");
 579     cell_index->SetNumberOfValues(grid->GetNumberOfCells());
 580     grid->GetCellData()->AddArray(cell_index);
 581   }
 582   EG_VTKDCC(vtkLongArray_t, cell_index, grid, "cell_index");
 583   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 584     cell_index->SetValue(id_cell, id_cell);
 585   }
 586 }
 587
 588 void EgVtkObject::UpdateNodeIndex(vtkUnstructuredGrid *grid)
 589 {
 590   if (!grid->GetPointData()->GetArray("node_index")) {
 591     EG_VTKSP(vtkLongArray_t, node_index);
 592     node_index->SetName("node_index");
 593     node_index->SetNumberOfValues(grid->GetNumberOfPoints());
 594     grid->GetPointData()->AddArray(node_index);
 595   }
 596   EG_VTKDCN(vtkLongArray_t, node_index, grid, "node_index");
 597   for (vtkIdType pointId = 0; pointId < grid->GetNumberOfPoints(); ++pointId) {
 598     node_index->SetValue(pointId, pointId);
 599   }
 600 }
 601
 602 void EgVtkObject::addToPolyData
 603 (
 604   QVector<vtkIdType>  &cells,
 605   vtkPolyData         *pdata,
 606   vtkUnstructuredGrid *grid
 607 )
 608 {
 609   UpdateCellIndex(grid);
 610   UpdateNodeIndex(grid);
 611   QVector<vtkIdType> nodes;
 612   QVector<int>       _nodes;
 613   getNodesFromCells(cells, nodes, grid);
 614   createNodeMapping(nodes, _nodes, grid);
 615   EG_VTKSP(vtkDoubleArray, pcoords);
 616   pcoords->SetNumberOfComponents(3);
 617   pcoords->SetNumberOfTuples(nodes.size());
 618   EG_VTKSP(vtkPoints, points);
 619   points->SetData(pcoords);
 620   pdata->SetPoints(points);
 621   pdata->Allocate(cells.size());
 622   if (!pdata->GetCellData()->GetArray("cell_index")) {
 623     EG_VTKSP(vtkLongArray_t, cell_index);
 624     cell_index->SetName("cell_index");
 625     //cell_index->SetNumberOfValues(cells.size());
 626     pdata->GetCellData()->AddArray(cell_index);
 627   }
 628   if (!pdata->GetPointData()->GetArray("node_index")) {
 629     EG_VTKSP(vtkLongArray_t, node_index);
 630     node_index->SetName("node_index");
 631     //node_index->SetNumberOfValues(nodes.size());
 632     pdata->GetPointData()->AddArray(node_index);
 633   }
 634   EG_VTKDCC(vtkLongArray_t, pd_cell_index, pdata, "cell_index");
 635   EG_VTKDCN(vtkLongArray_t, pd_node_index, pdata, "node_index");
 636   pd_cell_index->SetNumberOfValues(cells.size());
 637   pd_node_index->SetNumberOfValues(nodes.size());
 638   for (int i_cell = 0; i_cell < cells.size(); ++i_cell) {
 639     vtkIdType id_cell = cells[i_cell];
 640     vtkIdType cellType = grid->GetCellType(id_cell);
 641     if ((cellType != VTK_TRIANGLE) && (cellType != VTK_QUAD)) {
 642       EG_ERR_RETURN("unsupported cell type for this operation");
 643     }
 644     vtkIdType Npts, *pts;
 645     grid->GetCellPoints(id_cell, Npts, pts);
 646     vtkIdType *new_pts = new vtkIdType[Npts];
 647     for (int i = 0; i < Npts; ++i) {
 648       new_pts[i] = _nodes[pts[i]];
 649     }
 650     vtkIdType newCellId = pdata->InsertNextCell(cellType, Npts, new_pts);
 651     pd_cell_index->SetValue(newCellId, id_cell);
 652     delete [] new_pts;
 653   }
 654   for (int i_node = 0; i_node < nodes.size(); ++i_node) {
 655     vec3_t x;
 656     grid->GetPoints()->GetPoint(nodes[i_node], x.data());
 657     pdata->GetPoints()->SetPoint(i_node, x.data());
 658     pd_node_index->SetValue(i_node, nodes[i_node]);
 659   }
 660 }
 661
 662 #define EGVTKOBJECT_COPYCELLDATA(FIELD,TYPE) \
 663 { \
 664 if (old_grid->GetCellData()->GetArray(FIELD)) { \
 665 EG_VTKDCC(TYPE, var1, old_grid, FIELD); \
 666 EG_VTKDCC(TYPE, var2, new_grid, FIELD); \
 667 var2->SetValue(newId, var1->GetValue(oldId)); \
 668 } \
 669 }
 670
 671 void EgVtkObject::copyCellData
 672 (
 673   vtkUnstructuredGrid *old_grid,
 674   vtkIdType            oldId,
 675   vtkUnstructuredGrid *new_grid,
 676   vtkIdType            newId
 677 )
 678 {
 679   EGVTKOBJECT_COPYCELLDATA("vtk_type",    vtkIntArray);
 680   EGVTKOBJECT_COPYCELLDATA("cell_code",   vtkIntArray);
 681   EGVTKOBJECT_COPYCELLDATA("cell_orgdir", vtkIntArray);
 682   EGVTKOBJECT_COPYCELLDATA("cell_curdir", vtkIntArray);
 683   EGVTKOBJECT_COPYCELLDATA("cell_voldir", vtkIntArray);
 684   EGVTKOBJECT_COPYCELLDATA("cell_index",  vtkLongArray_t);
 685 }
 686
 687 #define EGVTKOBJECT_COPYNODEDATA(FIELD,TYPE) \
 688 { \
 689 if (old_grid->GetPointData()->GetArray(FIELD)) { \
 690 EG_VTKDCN(TYPE, var1, old_grid, FIELD); \
 691 EG_VTKDCN(TYPE, var2, new_grid, FIELD); \
 692 var2->SetValue(newId, var1->GetValue(oldId)); \
 693 } \
 694 }
 695
 696 void EgVtkObject::copyNodeData
 697 (
 698   vtkUnstructuredGrid *old_grid,
 699   vtkIdType            oldId,
 700   vtkUnstructuredGrid *new_grid,
 701   vtkIdType            newId
 702 )
 703 {
 704   EGVTKOBJECT_COPYNODEDATA("node_status", vtkIntArray);
 705   EGVTKOBJECT_COPYNODEDATA("node_layer", vtkIntArray);
 706   EGVTKOBJECT_COPYNODEDATA("node_index", vtkLongArray_t);
 707   EGVTKOBJECT_COPYNODEDATA("node_specified_density", vtkIntArray);
 708   EGVTKOBJECT_COPYNODEDATA("node_meshdensity_desired",  vtkDoubleArray);
 709   EGVTKOBJECT_COPYNODEDATA("node_meshdensity_current",  vtkDoubleArray);
 710   EGVTKOBJECT_COPYNODEDATA("node_type",  vtkCharArray);
 711 }
 712
 713 #define EGVTKOBJECT_CREATECELLFIELD(FIELD,TYPE,OW) \
 714 if (!grid->GetCellData()->GetArray(FIELD)) { \
 715   EG_VTKSP(TYPE, var); \
 716   var->SetName(FIELD); \
 717   var->SetNumberOfValues(Ncells); \
 718   grid->GetCellData()->AddArray(var); \
 719   for (int i = 0; i < grid->GetNumberOfCells(); ++i) { \
 720     var->SetValue(i,0); \
 721   } \
 722 } else if (OW) { \
 723   EG_VTKDCC(TYPE, var, grid, FIELD); \
 724   var->SetNumberOfValues(Ncells); \
 725   for (int i = 0; i < grid->GetNumberOfCells(); ++i) { \
 726     var->SetValue(i,0); \
 727   } \
 728 }
 729
 730 #define EGVTKOBJECT_CREATENODEFIELD(FIELD,TYPE,OW) \
 731 if (!grid->GetPointData()->GetArray(FIELD)) { \
 732   EG_VTKSP(TYPE, var); \
 733   var->SetName(FIELD); \
 734   var->SetNumberOfValues(Nnodes); \
 735   grid->GetPointData()->AddArray(var); \
 736   for (int i = 0; i < grid->GetNumberOfPoints(); ++i) { \
 737     var->SetValue(i,0); \
 738   } \
 739 } else if (OW) { \
 740   EG_VTKDCN(TYPE, var, grid, FIELD); \
 741   var->SetNumberOfValues(Nnodes); \
 742   for (int i = 0; i < grid->GetNumberOfPoints(); ++i) { \
 743     var->SetValue(i,0); \
 744   } \
 745 }
 746
 747 void EgVtkObject::createBasicFields(vtkUnstructuredGrid *grid, vtkIdType Ncells, vtkIdType Nnodes, bool overwrite)
 748 {
 749   createBasicNodeFields(grid, Nnodes, overwrite);
 750   createBasicCellFields(grid, Ncells, overwrite);
 751 }
 752
 753 void EgVtkObject::createBasicCellFields(vtkUnstructuredGrid *grid, vtkIdType Ncells, bool overwrite)
 754 {
 755   EGVTKOBJECT_CREATECELLFIELD("vtk_type" ,   vtkIntArray, overwrite);
 756   EGVTKOBJECT_CREATECELLFIELD("cell_code",   vtkIntArray, overwrite);
 757   EGVTKOBJECT_CREATECELLFIELD("cell_index",  vtkLongArray_t, overwrite);
 758   EGVTKOBJECT_CREATECELLFIELD("cell_orgdir", vtkIntArray, overwrite); // original orientation
 759   EGVTKOBJECT_CREATECELLFIELD("cell_curdir", vtkIntArray, overwrite); // current orientation
 760   EGVTKOBJECT_CREATECELLFIELD("cell_voldir", vtkIntArray, overwrite); // volume orientation -- only valid for a single (i.e. the current) volume
 761   EGVTKOBJECT_CREATECELLFIELD("cell_VA",     vtkDoubleArray, overwrite);
 762 }
 763
 764 void EgVtkObject::createBasicNodeFields(vtkUnstructuredGrid *grid, vtkIdType Nnodes, bool overwrite)
 765 {
 766   EGVTKOBJECT_CREATENODEFIELD("node_status",               vtkIntArray,    overwrite);
 767   EGVTKOBJECT_CREATENODEFIELD("node_layer",                vtkIntArray,    overwrite);
 768   EGVTKOBJECT_CREATENODEFIELD("node_index",                vtkLongArray_t, overwrite);
 769   EGVTKOBJECT_CREATENODEFIELD("node_specified_density",    vtkIntArray,    overwrite); //density index from table
 770   EGVTKOBJECT_CREATENODEFIELD("node_meshdensity_desired",  vtkDoubleArray, overwrite); //what we want
 771   EGVTKOBJECT_CREATENODEFIELD("node_meshdensity_current",  vtkDoubleArray, overwrite); //what we have
 772   EGVTKOBJECT_CREATENODEFIELD("node_type",                 vtkCharArray, overwrite);   //node type
 773 }
 774
 775 void EgVtkObject::allocateGrid(vtkUnstructuredGrid *grid, vtkIdType Ncells, vtkIdType Nnodes, bool create_fields)
 776 {
 777   EG_VTKSP(vtkPoints,points);
 778   points->SetNumberOfPoints(Nnodes);
 779   grid->SetPoints(points);
 780   grid->Allocate(Ncells,max(vtkIdType(1),Ncells/10));
 781   if (create_fields) {
 782     createBasicFields(grid, Ncells, Nnodes, true);
 783   }
 784 }
 785
 786 vec3_t EgVtkObject::cellCentre(vtkUnstructuredGrid *grid, vtkIdType id_cell)
 787 {
 788   vec3_t x,xc(0,0,0);
 789   vtkIdType *pts, N_pts;
 790   grid->GetCellPoints(id_cell, N_pts, pts);
 791   double f = 1.0/N_pts;
 792   for (int i_pts = 0; i_pts < N_pts; ++i_pts) {
 793     grid->GetPoint(pts[i_pts], x.data());
 794     xc += f*x;
 795   }
 796   return xc;
 797 }
 798
 799 void EgVtkObject::getRestCells(vtkUnstructuredGrid      *grid,
 800                                const QVector<vtkIdType> &cells,
 801                                QVector<vtkIdType>       &rest_cells)
 802 {
 803   QVector<bool> is_in_cells(grid->GetNumberOfCells(), false);
 804   foreach (vtkIdType id_cell, cells) {
 805     is_in_cells[id_cell] = true;
 806   }
 807   rest_cells.resize(grid->GetNumberOfCells() - cells.size());
 808   int i_rest_cells = 0;
 809   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 810     if (!is_in_cells[id_cell]) {
 811       rest_cells[i_rest_cells] = id_cell;
 812       ++i_rest_cells;
 813     }
 814   }
 815 }
 816
 817 void EgVtkObject::makeCopy(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst)
 818 {
 819   allocateGrid(dst, src->GetNumberOfCells(), src->GetNumberOfPoints());
 820   for (vtkIdType id_node = 0; id_node < src->GetNumberOfPoints(); ++id_node) {
 821     vec3_t x;
 822     src->GetPoints()->GetPoint(id_node, x.data());
 823     dst->GetPoints()->SetPoint(id_node, x.data());
 824     copyNodeData(src, id_node, dst, id_node);
 825   }
 826   for (vtkIdType id_cell = 0; id_cell < src->GetNumberOfCells(); ++id_cell) {
 827     vtkIdType N_pts, *pts;
 828     vtkIdType type_cell = src->GetCellType(id_cell);
 829     src->GetCellPoints(id_cell, N_pts, pts);
 830     vtkIdType id_new_cell = dst->InsertNextCell(type_cell, N_pts, pts);
 831     copyCellData(src, id_cell, dst, id_new_cell);
 832   }
 833 }
 834
 835 void EgVtkObject::makeCopyNoAlloc(vtkUnstructuredGrid *src, vtkUnstructuredGrid *dst)
 836 {
 837   for (vtkIdType id_node = 0; id_node < src->GetNumberOfPoints(); ++id_node) {
 838     vec3_t x;
 839     src->GetPoints()->GetPoint(id_node, x.data());
 840     dst->GetPoints()->SetPoint(id_node, x.data());
 841     copyNodeData(src, id_node, dst, id_node);
 842   }
 843   for (vtkIdType id_cell = 0; id_cell < src->GetNumberOfCells(); ++id_cell) {
 844     vtkIdType N_pts, *pts;
 845     vtkIdType type_cell = src->GetCellType(id_cell);
 846     src->GetCellPoints(id_cell, N_pts, pts);
 847     vtkIdType id_new_cell = dst->InsertNextCell(type_cell, N_pts, pts);
 848     copyCellData(src, id_cell, dst, id_new_cell);
 849   }
 850 }
 851
 852 void EgVtkObject::reorientateFace(vtkUnstructuredGrid *grid, vtkIdType id_face)
 853 {
 854   EG_VTKDCC(vtkIntArray, cell_curdir, grid, "cell_curdir");
 855   vtkIdType N_pts, *pts;
 856   grid->GetCellPoints(id_face, N_pts, pts);
 857   vtkIdType new_pts[N_pts];
 858   for (int i = 0; i < N_pts; ++i) {
 859     new_pts[i] = pts[N_pts - i - 1];
 860   }
 861   if (cell_curdir->GetValue(id_face) == 0) {
 862     cell_curdir->SetValue(id_face, 1);
 863   } else {
 864     cell_curdir->SetValue(id_face, 0);
 865   }
 866   grid->ReplaceCell(id_face, N_pts, new_pts);
 867 }
 868
 869 void EgVtkObject::resetOrientation(vtkUnstructuredGrid *grid)
 870 {
 871   EG_VTKDCC(vtkIntArray, cell_orgdir, grid, "cell_orgdir");
 872   EG_VTKDCC(vtkIntArray, cell_curdir, grid, "cell_curdir");
 873   EG_VTKDCC(vtkIntArray, cell_voldir, grid, "cell_voldir");
 874   QVector<vtkIdType> faces;
 875   getAllSurfaceCells(faces, grid);
 876   foreach (vtkIdType id_face, faces) {
 877     if (cell_curdir->GetValue(id_face) != cell_orgdir->GetValue(id_face)) {
 878       reorientateFace(grid, id_face);
 879       cell_curdir->SetValue(id_face, cell_orgdir->GetValue(id_face));
 880     }
 881     cell_voldir->SetValue(id_face, 0);
 882   }
 883 }
 884
 885 vtkIdType EgVtkObject::findVolumeCell(vtkUnstructuredGrid *grid, vtkIdType id_surf, g2l_t _nodes, l2g_t cells, g2l_t _cells, l2l_t n2c)
 886 {
 887   vtkIdType N_pts, *pts;
 888   if (_cells.size()) N_pts = N_pts; // dummy statement to get rid of compiler warning ...
 889   grid->GetCellPoints(id_surf, N_pts, pts);
 890   QVector<QSet<int> > inters(N_pts-1);
 891   qcontIntersection(n2c[_nodes[pts[0]]], n2c[_nodes[pts[1]]], inters[0]);
 892   int i_pts = 2;
 893   while (i_pts < N_pts) {
 894     qcontIntersection(inters[i_pts-2], n2c[_nodes[pts[i_pts]]], inters[i_pts-1]);
 895     ++i_pts;
 896   }
 897   if (inters[N_pts-2].size() == 0) {
 898     return -1;
 899   } else if (inters[N_pts-2].size() > 2) {
 900     EG_BUG;
 901   }
 902   vtkIdType id_vol = -1;
 903   foreach (int i_cells, inters[N_pts-2]) {
 904     if (cells[i_cells] != id_surf) {
 905       id_vol = cells[i_cells];
 906     }
 907   }
 908   return id_vol;
 909 }
 910
 911 ///\todo Why not simply use m_BoundaryCodes = bcs; ?
 912 //WARNING: Do never call this->setBoundaryCodes(this->m_BoundaryCodes) as this will empty m_BoundaryCodes!!!
 913 void EgVtkObject::setBoundaryCodes(const QSet<int> &bcs)
 914 {
 915   m_BoundaryCodes = bcs;
 916 /*  m_BoundaryCodes.clear();
 917   int bc;
 918   foreach(bc, bcs) {
 919     m_BoundaryCodes.insert(bc);
 920   }*/
 921 }
 922
 923 QSet<int> EgVtkObject::getBoundaryCodes()
 924 {
 925   return m_BoundaryCodes;
 926 }
 927
 928 void EgVtkObject::createIndices(vtkUnstructuredGrid *grid)
 929 {
 930   if (!grid->GetCellData()->GetArray("cell_index")) {
 931     EG_VTKSP(vtkLongArray_t, var);
 932     var->SetName("cell_index");
 933     var->SetNumberOfValues(grid->GetNumberOfCells());
 934     grid->GetCellData()->AddArray(var);
 935   } else {
 936     EG_VTKDCC(vtkLongArray_t, var, grid, "cell_index");
 937     var->SetNumberOfValues(grid->GetNumberOfCells());
 938   }
 939   EG_VTKDCC(vtkLongArray_t, cell_index, grid, "cell_index");
 940   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 941     cell_index->SetValue(id_cell, id_cell);
 942   }
 943
 944   if (!grid->GetCellData()->GetArray("vtk_type")) {
 945     EG_VTKSP(vtkIntArray, var);
 946     var->SetName("vtk_type");
 947     var->SetNumberOfValues(grid->GetNumberOfCells());
 948     grid->GetCellData()->AddArray(var);
 949   } else {
 950     EG_VTKDCC(vtkIntArray, var, grid, "vtk_type");
 951     var->SetNumberOfValues(grid->GetNumberOfCells());
 952   }
 953   EG_VTKDCC(vtkIntArray, vtk_type, grid, "vtk_type");
 954   for (vtkIdType id_cell = 0; id_cell < grid->GetNumberOfCells(); ++id_cell) {
 955     vtk_type->SetValue(id_cell, grid->GetCellType(id_cell));
 956   }
 957
 958   if (!grid->GetCellData()->GetArray("node_index")) {
 959     EG_VTKSP(vtkLongArray_t, var);
 960     var->SetName("node_index");
 961     var->SetNumberOfValues(grid->GetNumberOfPoints());
 962     grid->GetPointData()->AddArray(var);
 963   } else {
 964     EG_VTKDCC(vtkLongArray_t, var, grid, "node_index");
 965     var->SetNumberOfValues(grid->GetNumberOfPoints());
 966   }
 967   EG_VTKDCN(vtkLongArray_t, node_index, grid, "node_index");
 968   for (vtkIdType id_node = 0; id_node < grid->GetNumberOfPoints(); ++id_node) {
 969     node_index->SetValue(id_node, id_node);
 970   }
 971 }
 972
 973 BoundaryCondition EgVtkObject::getBC(int bc)
 974 {
 975   return GuiMainWindow::pointer()->getBC(bc);
 976 }
 977
 978 int EgVtkObject::getSet(QString group, QString key, int value, int& variable)
 979 {
 980   QSettings *qset = GuiMainWindow::settings();
 981   QString typed_key = "int/" + key;
 982   if(group!=QObject::tr("General")) qset->beginGroup(group);
 983   //if key=value pair not found in settings file, write it
 984   if (!qset->contains(typed_key)) qset->setValue(typed_key,value);
 985   //read key value from settings file and assign it to variable
 986   variable = (qset->value(typed_key,variable)).toInt();
 987   if(group!=QObject::tr("General")) qset->endGroup();
 988   return(variable);
 989 }
 990
 991 double EgVtkObject::getSet(QString group, QString key, double value, double& variable)
 992 {
 993   QSettings *qset = GuiMainWindow::settings();
 994   QString typed_key = "double/" + key;
 995   if(group!=QObject::tr("General")) qset->beginGroup(group);
 996   //if key=value pair not found in settings file, write it
 997   if (!qset->contains(typed_key)) qset->setValue(typed_key,value);
 998   //read key value from settings file and assign it to variable
 999   variable = (qset->value(typed_key,variable)).toDouble();
1000   if(group!=QObject::tr("General")) qset->endGroup();
1001   return(variable);
1002 }
1003
1004 bool EgVtkObject::getSet(QString group, QString key, bool value, bool& variable)
1005 {
1006   QSettings *qset = GuiMainWindow::settings();
1007   QString typed_key = "bool/" + key;
1008   if(group!=QObject::tr("General")) qset->beginGroup(group);
1009   Qt::CheckState state = (Qt::CheckState) ( value ? 2 : 0 );
1010   //if key=value pair not found in settings file, write it
1011   if (!qset->contains(typed_key)) qset->setValue(typed_key,state);
1012   //read key value from settings file and assign it to variable
1013   variable = (qset->value(typed_key,variable)).toBool();
1014   if(group!=QObject::tr("General")) qset->endGroup();
1015   return(variable);
1016 }
1017
1018 QString EgVtkObject::getSet(QString group, QString key, QString value, QString& variable)
1019 {
1020   QSettings *qset = GuiMainWindow::settings();
1021   QString typed_key;
1022   typed_key = QObject::tr("QString/") + key;
1023   if (group != QObject::tr("General")) qset->beginGroup(group);
1024   //if key=value pair not found in settings file, write it
1025   if (!qset->contains(typed_key)) qset->setValue(typed_key, value);
1026   //read key value from settings file and assign it to variable
1027   variable = (qset->value(typed_key)).toString();
1028   if (group != QObject::tr("General")) qset->endGroup();
1029   return(variable);
1030 }
1031
1032 QString EgVtkObject::getSet(QString group, QString key, QString value, QString& variable, int type)
1033 {
1034   QSettings *qset = GuiMainWindow::settings();
1035   QString typed_key;
1036   if (type == 0) {
1037     typed_key = QObject::tr("QString/") + key;
1038   }
1039   else if (type == 1) {
1040     typed_key = QObject::tr("Filename/") + key;
1041   }
1042   else {
1043     typed_key = QObject::tr("Directory/") + key;
1044   }
1045   if (group != QObject::tr("General")) qset->beginGroup(group);
1046   //if key=value pair not found in settings file, write it
1047   if (!qset->contains(typed_key)) qset->setValue(typed_key, value);
1048   //read key value from settings file and assign it to variable
1049   variable = (qset->value(typed_key)).toString();
1050   if (group != QObject::tr("General")) qset->endGroup();
1051   return(variable);
1052 }
1053
1054 void EgVtkObject::writeGrid(vtkUnstructuredGrid *grid, QString name)
1055 {
1056   QVector<vtkIdType> cells;
1057   getAllCells(cells, grid);
1058   name = GuiMainWindow::pointer()->getCwd() + "/" + name + ".vtu";
1059   writeCells(grid, cells, name);
1060 //   qDebug()<<"Saved grid as "<<name;
1061 }
1062
1063 void EgVtkObject::getAllNodeDataNames(QVector<QString> &field_names, vtkUnstructuredGrid *grid)
1064 {
1065   int N = grid->GetPointData()->GetNumberOfArrays();
1066   field_names.resize(N);
1067   for (int i = 0; i < N; ++i) {
1068     field_names[i] = grid->GetPointData()->GetArrayName(i);
1069   }
1070 }
1071
1072 void EgVtkObject::getAllCellDataNames(QVector<QString> &field_names, vtkUnstructuredGrid *grid)
1073 {
1074   int N = grid->GetCellData()->GetNumberOfArrays();
1075   field_names.resize(N);
1076   for (int i = 0; i < N; ++i) {
1077     field_names[i] = grid->GetCellData()->GetArrayName(i);
1078   }
1079 }
1080
1081 QString EgVtkObject::stripFromExtension(QString file_name)
1082 {
1083   int i = file_name.size() - 1;
1084   while ((i > 0) && (file_name[i] != '.') && (file_name[i] != '/') && (file_name[i] != '\\')) {
1085     --i;
1086   }
1087   if (file_name[i] == '.') {
1088     return file_name.left(i);
1089   }
1090   return file_name;
1091 }
1092
1093 QString EgVtkObject::getExtension(QString file_name)
1094 {
1095   int i = file_name.size();
1096   while ((i > 0) && (file_name[i] != '.') && (file_name[i] != '/') && (file_name[i] != '\\')) {
1097     --i;
1098   }
1099   if (file_name[i] == '.') {
1100     return (file_name.right(file_name.size() - i - 1)).toLower();
1101   }
1102   return "";
1103 }
1104
1105 ///////////////////////////////////////////
1106
1107 void EgVtkObject::getFaceOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_face, QVector<vtkIdType> &ids)
1108 {
1109   vtkIdType type_cell = grid->GetCellType(id_cell);
1110   ids.clear();
1111   vtkIdType *pts, N_pts;
1112   grid->GetCellPoints(id_cell, N_pts, pts);
1113   if (type_cell == VTK_TETRA) {
1114     ids.resize(3);
1115     if      (i_face == 0) { ids[0] = pts[2]; ids[1] = pts[1]; ids[2] = pts[0]; }
1116     else if (i_face == 1) { ids[0] = pts[1]; ids[1] = pts[3]; ids[2] = pts[0]; }
1117     else if (i_face == 2) { ids[0] = pts[3]; ids[1] = pts[2]; ids[2] = pts[0]; }
1118     else if (i_face == 3) { ids[0] = pts[2]; ids[1] = pts[3]; ids[2] = pts[1]; }
1119   } else {
1120     EG_BUG; // not implemented
1121   }
1122 }
1123
1124 void EgVtkObject::getEdgeOfCell(vtkUnstructuredGrid *grid, vtkIdType id_cell, int i_edge, QVector<vtkIdType> &ids)
1125 {
1126   vtkIdType type_cell = grid->GetCellType(id_cell);
1127   ids.clear();
1128   vtkIdType *pts, N_pts;
1129   grid->GetCellPoints(id_cell, N_pts, pts);
1130   if (type_cell == VTK_TETRA) {
1131     ids.resize(2);
1132     if      (i_edge == 0) { ids[0] = pts[0]; ids[1] = pts[1]; }
1133     else if (i_edge == 1) { ids[0] = pts[0]; ids[1] = pts[2]; }
1134     else if (i_edge == 2) { ids[0] = pts[0]; ids[1] = pts[3]; }
1135     else if (i_edge == 3) { ids[0] = pts[1]; ids[1] = pts[2]; }
1136     else if (i_edge == 4) { ids[0] = pts[1]; ids[1] = pts[3]; }
1137     else if (i_edge == 5) { ids[0] = pts[2]; ids[1] = pts[3]; }
1138   } else {
1139     EG_BUG; // not implemented
1140   }
1141 }
1142
1143 bool EgVtkObject::saveGrid(vtkUnstructuredGrid* a_grid, QString file_name)
1144 {
1145   file_name += ".vtu";
1146   addVtkTypeInfo(a_grid);
1147   createIndices(a_grid);
1148   EG_VTKSP(vtkXMLUnstructuredGridWriter,vtu);
1149   vtu->SetFileName(qPrintable(file_name));
1150   vtu->SetDataModeToBinary();
1151   vtu->SetInput(a_grid);
1152   vtu->Write();
1153   if(vtu->GetErrorCode()) {
1154     return false;
1155   }
1156   else {
1157     return true;
1158   }
1159 }
1160
1161 void EgVtkObject::addVtkTypeInfo(vtkUnstructuredGrid* a_grid)
1162 {
1163   EG_VTKSP(vtkIntArray, vtk_type);
1164   vtk_type->SetName("vtk_type");
1165   vtk_type->SetNumberOfValues(a_grid->GetNumberOfCells());
1166   EG_VTKDCC(vtkDoubleArray, cell_VA, a_grid, "cell_VA");
1167   for (vtkIdType cellId = 0; cellId < a_grid->GetNumberOfCells(); ++cellId) {
1168     vtk_type->SetValue(cellId, a_grid->GetCellType(cellId));
1169     cell_VA->SetValue(cellId, GeometryTools::cellVA(a_grid, cellId, true));
1170   }
1171   a_grid->GetCellData()->AddArray(vtk_type);
1172 }
1173
1174 vtkIdType EgVtkObject::addGrid(vtkUnstructuredGrid *main_grid, vtkUnstructuredGrid *grid_to_add, vtkIdType offset)
1175 {
1176   for (vtkIdType id_node = 0; id_node < grid_to_add->GetNumberOfPoints(); ++id_node) {
1177     vec3_t x;
1178     grid_to_add->GetPoints()->GetPoint(id_node, x.data());
1179     main_grid->GetPoints()->SetPoint(offset + id_node, x.data());
1180     copyNodeData(grid_to_add, id_node, main_grid, offset + id_node);
1181   }
1182   for (vtkIdType id_cell = 0; id_cell < grid_to_add->GetNumberOfCells(); ++id_cell) {
1183     vtkIdType N_pts, *pts;
1184     vtkIdType type_cell = grid_to_add->GetCellType(id_cell);
1185     grid_to_add->GetCellPoints(id_cell, N_pts, pts);
1186     QVector <vtkIdType> new_pts(N_pts);
1187     for(int i=0;i<N_pts;i++) new_pts[i] = offset + pts[i];
1188     vtkIdType id_new_cell = main_grid->InsertNextCell(type_cell, N_pts, new_pts.data());
1189     copyCellData(grid_to_add, id_cell, main_grid, id_new_cell);
1190   }
1191   return( offset + grid_to_add->GetNumberOfPoints() );
1192 }
1193