src/VTKViewer/VTKViewer_GeometryFilter.cxx

   1 // Copyright (C) 2007-2013  CEA/DEN, EDF R&D, OPEN CASCADE
   2 //
   3 // Copyright (C) 2003-2007  OPEN CASCADE, EADS/CCR, LIP6, CEA/DEN,
   4 // CEDRAT, EDF R&D, LEG, PRINCIPIA R&D, BUREAU VERITAS
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 //  File   : VTKViewer_GeometryFilter.cxx
  24 //  Author : Michael ZORIN
  25 //  Module : SALOME
  26 //
  27 #include "VTKViewer_GeometryFilter.h"
  28 #include "VTKViewer_ConvexTool.h"
  29 #include "VTKViewer_ArcBuilder.h"
  30
  31 #include <vtkSmartPointer.h>
  32 #include <vtkCellArray.h>
  33 #include <vtkCellData.h>
  34 #include <vtkGenericCell.h>
  35 #include <vtkHexagonalPrism.h>
  36 #include <vtkHexahedron.h>
  37 #include <vtkInformation.h>
  38 #include <vtkInformationVector.h>
  39 #include <vtkMergePoints.h>
  40 #include <vtkObjectFactory.h>
  41 #include <vtkPointData.h>
  42 #include <vtkPolyData.h>
  43 #include <vtkPolygon.h>
  44 #include <vtkPyramid.h>
  45 #include <vtkStructuredGrid.h>
  46 #include <vtkTetra.h>
  47 #include <vtkUnsignedCharArray.h>
  48 #include <vtkUnstructuredGrid.h>
  49 #include <vtkVoxel.h>
  50 #include <vtkWedge.h>
  51
  52 #include <algorithm>
  53 #include <iterator>
  54 #include <vector>
  55 #include <map>
  56 #include <set>
  57
  58 #include "utilities.h"
  59
  60 #if defined __GNUC__
  61   #if __GNUC__ == 2
  62     #define __GNUC_2__
  63   #endif
  64 #endif
  65
  66 #define VTK_XVERSION (VTK_MAJOR_VERSION*10000+VTK_MINOR_VERSION*100+VTK_BUILD_VERSION)
  67
  68 //#define __MYDEBUG__
  69 //#define USE_ROBUST_TRIANGULATION
  70
  71 ///////////////////////////////////////////////////////////////////////////////////////////////
  72 // VSR 26/10/2012: fix of regression (issue 21924) - increased memory consumption
  73 // for displaying of 3d elements, introduced by fix for issue 20314.
  74 // ...
  75 // The macro SHOW_COINCIDING_3D_PAL20314, when defined, allows correct visualization of
  76 // coincident 3d elements but causes substantial increasing of memory consumption, as all 3d
  77 // elements are always shown, even if they are totally covered by surrounding faces.
  78 // If this macro is not defined (commented), the behaviour is defined by another macro -
  79 // SHOW_COINCIDING_3D_PAL21924, as follows:
  80 // - If SHOW_COINCIDING_3D_PAL21924 is defined, an alternative solution for computing
  81 //   visibility of 3d elements is used; this solution allows to fix problem with visibility
  82 //   of coinciding 3d elements in most cases (though some cases might not work), while not
  83 //   causing significant increasing of memory consumption.
  84 // - If SHOW_COINCIDING_3D_PAL21924 is not defined (commented), coinciding 3d elements are
  85 //   not shown at all (this corresponds to the state before issue 20314 fixing).
  86 ///////////////////////////////////////////////////////////////////////////////////////////////
  87 //#define SHOW_COINCIDING_3D_PAL20314
  88 #ifndef SHOW_COINCIDING_3D_PAL20314
  89 #define SHOW_COINCIDING_3D_PAL21924
  90 #endif
  91 ///////////////////////////////////////////////////////////////////////////////////////////////
  92
  93 vtkStandardNewMacro(VTKViewer_GeometryFilter);
  94
  95 VTKViewer_GeometryFilter
  96 ::VTKViewer_GeometryFilter():
  97   myShowInside(0),
  98   myStoreMapping(0),
  99   myAppendCoincident3D(0),
 100   myIsWireframeMode(0),
 101   myIsBuildArc(false),
 102   myMaxArcAngle(2)
 103 {}
 104
 105
 106 VTKViewer_GeometryFilter
 107 ::~VTKViewer_GeometryFilter()
 108 {}
 109
 110
 111 int
 112 VTKViewer_GeometryFilter
 113 ::RequestData(
 114   vtkInformation *request,
 115   vtkInformationVector **inputVector,
 116   vtkInformationVector *outputVector)
 117 {
 118   // get the info objects
 119   vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
 120   vtkInformation *outInfo = outputVector->GetInformationObject(0);
 121
 122   // get the input and ouptut
 123   vtkDataSet *input = vtkDataSet::SafeDownCast(
 124     inInfo->Get(vtkDataObject::DATA_OBJECT()));
 125   vtkPolyData *output = vtkPolyData::SafeDownCast(
 126     outInfo->Get(vtkDataObject::DATA_OBJECT()));
 127
 128   vtkIdType numCells=input->GetNumberOfCells();
 129
 130   if (numCells == 0)
 131     {
 132       return 0;
 133     }
 134
 135   if (input->GetDataObjectType() == VTK_UNSTRUCTURED_GRID){
 136     return this->UnstructuredGridExecute(input, output, outInfo);
 137   }else
 138     return Superclass::RequestData(request,inputVector,outputVector);
 139
 140   return 1;
 141 }
 142
 143
 144 int
 145 VTKViewer_GeometryFilter
 146 ::UnstructuredGridExecute(
 147                           vtkDataSet *dataSetInput,
 148                           vtkPolyData *output,
 149                           vtkInformation *outInfo)
 150 {
 151
 152   vtkUnstructuredGrid *input= (vtkUnstructuredGrid *)dataSetInput;
 153   vtkCellArray *Connectivity = input->GetCells();
 154   // Check input
 155   if ( Connectivity == NULL )
 156     {
 157     vtkDebugMacro(<<"Nothing to extract");
 158     return 0;
 159     }
 160
 161   vtkIdType cellId;
 162   int i;
 163   int allVisible;
 164   vtkIdType npts = 0;
 165   vtkIdType *pts = 0;
 166   vtkPoints *p = input->GetPoints();
 167   vtkIdType numCells=input->GetNumberOfCells();
 168   vtkPointData *pd = input->GetPointData();
 169   vtkCellData *cd = input->GetCellData();
 170   vtkPointData *outputPD = output->GetPointData();
 171
 172   VTKViewer_OrderedTriangulator anOrderedTriangulator;
 173   VTKViewer_DelaunayTriangulator aDelaunayTriangulator;
 174
 175   vtkCellData *outputCD = output->GetCellData();
 176   vtkGenericCell *cell = vtkGenericCell::New();
 177
 178   vtkIdList *cellIds = vtkIdList::New();
 179   vtkIdList *faceIds = vtkIdList::New();
 180   vtkIdList *cellIdsTmp = vtkIdList::New();
 181   vtkIdList *faceIdsTmp = vtkIdList::New();
 182
 183   char *cellVis;
 184   vtkIdType newCellId;
 185   int faceId, *faceVerts, numFacePts;
 186   double *x;
 187   vtkIdType PixelConvert[4];
 188   // Change the type from int to vtkIdType in order to avoid compilation errors while using VTK
 189   // from ParaView-3.4.0 compiled on 64-bit Debian platform with VTK_USE_64BIT_IDS = ON
 190   vtkIdType aNewPts[VTK_CELL_SIZE];
 191   // ghost cell stuff
 192   unsigned char  updateLevel = (unsigned char)(GetUpdateGhostLevel());
 193   unsigned char  *cellGhostLevels = 0;
 194
 195   PixelConvert[0] = 0;
 196   PixelConvert[1] = 1;
 197   PixelConvert[2] = 3;
 198   PixelConvert[3] = 2;
 199
 200   vtkDebugMacro(<<"Executing geometry filter for unstructured grid input");
 201
 202   vtkDataArray* temp = 0;
 203   if (cd)
 204     {
 205     temp = cd->GetArray("vtkGhostLevels");
 206     }
 207   if ( (!temp) || (temp->GetDataType() != VTK_UNSIGNED_CHAR)
 208     || (temp->GetNumberOfComponents() != 1))
 209     {
 210     vtkDebugMacro("No appropriate ghost levels field available.");
 211     }
 212   else
 213     {
 214     cellGhostLevels = ((vtkUnsignedCharArray*)temp)->GetPointer(0);
 215     }
 216
 217   // Determine nature of what we have to do
 218   if ( (!this->CellClipping) && (!this->PointClipping) &&
 219        (!this->ExtentClipping) )
 220     {
 221     allVisible = 1;
 222     cellVis = NULL;
 223     }
 224   else
 225     {
 226     allVisible = 0;
 227     cellVis = new char[numCells];
 228     }
 229
 230   // Issue 0020115: [CEA 308] Quadratic elements visualization
 231   // Fix of remark described in note 0005222 - SIGSEGV
 232   vtkPoints* outputPoints = vtkPoints::New();
 233   outputPoints->DeepCopy(input->GetPoints());
 234   output->SetPoints(outputPoints);
 235   outputPoints->Delete();
 236
 237   outputPD->PassData(pd);
 238
 239   outputCD->CopyAllocate(cd,numCells,numCells/2);
 240
 241   output->Allocate(numCells/4+1,numCells);
 242
 243   // Loop over the cells determining what's visible
 244   if (!allVisible)
 245     {
 246     for (cellId=0, Connectivity->InitTraversal();
 247          Connectivity->GetNextCell(npts,pts);
 248          cellId++)
 249       {
 250       cellVis[cellId] = 1;
 251       if ( ( this->CellClipping && cellId < this->CellMinimum ) ||
 252            cellId > this->CellMaximum )
 253         {
 254         cellVis[cellId] = 0;
 255         }
 256       else
 257         {
 258         for (i=0; i < npts; i++)
 259           {
 260           x = p->GetPoint(pts[i]);
 261           if ( ( ( ( this->PointClipping && (pts[i] < this->PointMinimum ) ) ||
 262                                              pts[i] > this->PointMaximum) ) ||
 263                ( this->ExtentClipping &&
 264                 ( x[0] < this->Extent[0] || x[0] > this->Extent[1] ||
 265                   x[1] < this->Extent[2] || x[1] > this->Extent[3] ||
 266                   x[2] < this->Extent[4] || x[2] > this->Extent[5] )) )
 267             {
 268             cellVis[cellId] = 0;
 269             break;
 270             }//point/extent clipping
 271           }//for each point
 272         }//if point clipping needs checking
 273       }//for all cells
 274     }//if not all visible
 275
 276   // Loop over all cells now that visibility is known
 277   // (Have to compute visibility first for 3D cell boundarys)
 278   int progressInterval = numCells/20 + 1;
 279   TMapOfVectorId aDimension2VTK2ObjIds;
 280   if(myStoreMapping){
 281     myVTK2ObjIds.clear();
 282     myVTK2ObjIds.reserve(numCells);
 283   }
 284   for (cellId=0, Connectivity->InitTraversal();
 285        Connectivity->GetNextCell(npts,pts);
 286        cellId++)
 287     {
 288     //Progress and abort method support
 289     if ( !(cellId % progressInterval) )
 290       {
 291       vtkDebugMacro(<<"Process cell #" << cellId);
 292       this->UpdateProgress ((float)cellId/numCells);
 293       }
 294
 295     // Handle ghost cells here.  Another option was used cellVis array.
 296     if (cellGhostLevels && cellGhostLevels[cellId] > updateLevel)
 297       { // Do not create surfaces in outer ghost cells.
 298       continue;
 299       }
 300
 301     if (allVisible || cellVis[cellId])  //now if visible extract geometry
 302       {
 303       //special code for nonlinear cells - rarely occurs, so right now it
 304       //is slow.
 305       vtkIdType aCellType = input->GetCellType(cellId);
 306       switch (aCellType)
 307         {
 308         case VTK_EMPTY_CELL:
 309           break;
 310
 311         case VTK_VERTEX:
 312         case VTK_POLY_VERTEX:
 313           newCellId = output->InsertNextCell(aCellType,npts,pts);
 314           if(myStoreMapping){
 315             InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 316           }
 317           outputCD->CopyData(cd,cellId,newCellId);
 318           break;
 319
 320         case VTK_LINE:
 321         case VTK_POLY_LINE:
 322           newCellId = output->InsertNextCell(aCellType,npts,pts);
 323           if(myStoreMapping)
 324             InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 325           outputCD->CopyData(cd,cellId,newCellId);
 326           break;
 327
 328         case VTK_TRIANGLE:
 329         case VTK_QUAD:
 330         case VTK_POLYGON:
 331           newCellId = output->InsertNextCell(aCellType,npts,pts);
 332           if(myStoreMapping)
 333             InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 334           outputCD->CopyData(cd,cellId,newCellId);
 335           break;
 336
 337         case VTK_TRIANGLE_STRIP:
 338           newCellId = output->InsertNextCell(aCellType,npts,pts);
 339           if(myStoreMapping)
 340             InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 341           outputCD->CopyData(cd,cellId,newCellId);
 342           break;
 343
 344         case VTK_PIXEL:
 345           newCellId = output->InsertNextCell(aCellType,npts,pts);
 346           if(myStoreMapping)
 347             InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 348           outputCD->CopyData(cd,cellId,newCellId);
 349           break;
 350
 351         case VTK_CONVEX_POINT_SET: {
 352           bool anIsOk = anOrderedTriangulator.Execute(input,
 353                                                       cd,
 354                                                       cellId,
 355                                                       myShowInside,
 356                                                       allVisible,
 357                                                       GetAppendCoincident3D(),
 358                                                       cellVis,
 359                                                       output,
 360                                                       outputCD,
 361                                                       myStoreMapping,
 362                                                       myVTK2ObjIds,
 363                                                       aDimension2VTK2ObjIds,
 364                                                       true);
 365           if(!anIsOk)
 366             aDelaunayTriangulator.Execute(input,
 367                                           cd,
 368                                           cellId,
 369                                           myShowInside,
 370                                           allVisible,
 371                                           GetAppendCoincident3D(),
 372                                           cellVis,
 373                                           output,
 374                                           outputCD,
 375                                           myStoreMapping,
 376                                           myVTK2ObjIds,
 377                                           aDimension2VTK2ObjIds,
 378                                           false);
 379
 380           break;
 381         }
 382         case VTK_TETRA: {
 383 #ifdef SHOW_COINCIDING_3D_PAL21924
 384           faceIdsTmp->Reset();
 385           for (int ai=0; ai<npts; ai++)
 386             faceIdsTmp->InsertNextId(pts[ai]);
 387           input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 388 #endif
 389           for (faceId = 0; faceId < 4; faceId++)
 390             {
 391             faceIds->Reset();
 392             faceVerts = vtkTetra::GetFaceArray(faceId);
 393             faceIds->InsertNextId(pts[faceVerts[0]]);
 394             faceIds->InsertNextId(pts[faceVerts[1]]);
 395             faceIds->InsertNextId(pts[faceVerts[2]]);
 396             aCellType = VTK_TRIANGLE;
 397             numFacePts = 3;
 398             input->GetCellNeighbors(cellId, faceIds, cellIds);
 399 #ifdef SHOW_COINCIDING_3D_PAL21924
 400             int nbNeighbors = 0;
 401             for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 402               if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 403             }
 404             bool process = nbNeighbors <= 0;
 405 #else
 406             bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 407 #endif
 408             if ( process || myShowInside ||
 409                  (!allVisible && !cellVis[cellIds->GetId(0)]) )
 410               {
 411               for ( i=0; i < numFacePts; i++)
 412                 aNewPts[i] = pts[faceVerts[i]];
 413               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 414               if(myStoreMapping)
 415                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 416               outputCD->CopyData(cd,cellId,newCellId);
 417               }
 418             }
 419           break;
 420         }
 421         case VTK_VOXEL: {
 422 #ifdef SHOW_COINCIDING_3D_PAL21924
 423           faceIdsTmp->Reset();
 424           for (int ai=0; ai<npts; ai++)
 425             faceIdsTmp->InsertNextId(pts[ai]);
 426           input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 427 #endif
 428           for (faceId = 0; faceId < 6; faceId++)
 429             {
 430             faceIds->Reset();
 431             faceVerts = vtkVoxel::GetFaceArray(faceId);
 432             faceIds->InsertNextId(pts[faceVerts[0]]);
 433             faceIds->InsertNextId(pts[faceVerts[1]]);
 434             faceIds->InsertNextId(pts[faceVerts[2]]);
 435             faceIds->InsertNextId(pts[faceVerts[3]]);
 436             aCellType = VTK_QUAD;
 437             numFacePts = 4;
 438             input->GetCellNeighbors(cellId, faceIds, cellIds);
 439 #ifdef SHOW_COINCIDING_3D_PAL21924
 440             int nbNeighbors = 0;
 441             for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 442               if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 443             }
 444             bool process = nbNeighbors <= 0;
 445 #else
 446             bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 447 #endif
 448             if ( process || myShowInside ||
 449                  (!allVisible && !cellVis[cellIds->GetId(0)]) )
 450               {
 451               for ( i=0; i < numFacePts; i++)
 452                 aNewPts[i] = pts[faceVerts[PixelConvert[i]]];
 453               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 454               if(myStoreMapping)
 455                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 456               outputCD->CopyData(cd,cellId,newCellId);
 457               }
 458             }
 459           break;
 460         }
 461         case VTK_HEXAHEDRON: {
 462 #ifdef SHOW_COINCIDING_3D_PAL21924
 463           faceIdsTmp->Reset();
 464           for (int ai=0; ai<npts; ai++)
 465             faceIdsTmp->InsertNextId(pts[ai]);
 466           input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 467 #endif
 468           for (faceId = 0; faceId < 6; faceId++)
 469             {
 470             faceIds->Reset();
 471             faceVerts = vtkHexahedron::GetFaceArray(faceId);
 472             faceIds->InsertNextId(pts[faceVerts[0]]);
 473             faceIds->InsertNextId(pts[faceVerts[1]]);
 474             faceIds->InsertNextId(pts[faceVerts[2]]);
 475             faceIds->InsertNextId(pts[faceVerts[3]]);
 476             aCellType = VTK_QUAD;
 477             numFacePts = 4;
 478             input->GetCellNeighbors(cellId, faceIds, cellIds);
 479 #ifdef SHOW_COINCIDING_3D_PAL21924
 480             int nbNeighbors = 0;
 481             for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 482               if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 483             }
 484             bool process = nbNeighbors <= 0;
 485 #else
 486             bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 487 #endif
 488             if ( process || myShowInside ||
 489                  (!allVisible && !cellVis[cellIds->GetId(0)]) )
 490               {
 491               for ( i=0; i < numFacePts; i++)
 492                 aNewPts[i] = pts[faceVerts[i]];
 493               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 494               if(myStoreMapping)
 495                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 496               outputCD->CopyData(cd,cellId,newCellId);
 497               }
 498             }
 499           break;
 500         }
 501         case VTK_WEDGE: {
 502 #ifdef SHOW_COINCIDING_3D_PAL21924
 503           faceIdsTmp->Reset();
 504           for (int ai=0; ai<npts; ai++)
 505             faceIdsTmp->InsertNextId(pts[ai]);
 506           input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 507 #endif
 508           for (faceId = 0; faceId < 5; faceId++)
 509             {
 510             faceIds->Reset();
 511             faceVerts = vtkWedge::GetFaceArray(faceId);
 512             faceIds->InsertNextId(pts[faceVerts[0]]);
 513             faceIds->InsertNextId(pts[faceVerts[1]]);
 514             faceIds->InsertNextId(pts[faceVerts[2]]);
 515             aCellType = VTK_TRIANGLE;
 516             numFacePts = 3;
 517             if (faceVerts[3] >= 0)
 518               {
 519               faceIds->InsertNextId(pts[faceVerts[3]]);
 520               aCellType = VTK_QUAD;
 521               numFacePts = 4;
 522               }
 523
 524             input->GetCellNeighbors(cellId, faceIds, cellIds);
 525 #ifdef SHOW_COINCIDING_3D_PAL21924
 526             int nbNeighbors = 0;
 527             for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 528               if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 529             }
 530             bool process = nbNeighbors <= 0;
 531 #else
 532             bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 533 #endif
 534             if ( process || myShowInside ||
 535                  (!allVisible && !cellVis[cellIds->GetId(0)]) )
 536               {
 537               for ( i=0; i < numFacePts; i++)
 538                 aNewPts[i] = pts[faceVerts[i]];
 539               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 540               if(myStoreMapping)
 541                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 542               outputCD->CopyData(cd,cellId,newCellId);
 543               }
 544             }
 545           break;
 546         }
 547         case VTK_HEXAGONAL_PRISM: {
 548 #ifdef SHOW_COINCIDING_3D_PAL21924
 549           faceIdsTmp->Reset();
 550           for (int ai=0; ai<npts; ai++)
 551             faceIdsTmp->InsertNextId(pts[ai]);
 552           input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 553 #endif
 554           for (faceId = 0; faceId < 8; faceId++)
 555           {
 556             faceVerts = vtkHexagonalPrism::GetFaceArray(faceId);
 557             faceIds->Reset();
 558             faceIds->InsertNextId(pts[faceVerts[0]]);
 559             faceIds->InsertNextId(pts[faceVerts[1]]);
 560             faceIds->InsertNextId(pts[faceVerts[2]]);
 561             faceIds->InsertNextId(pts[faceVerts[3]]);
 562             aCellType = VTK_QUAD;
 563             numFacePts = 4;
 564             if (faceVerts[5] >= 0)
 565             {
 566               faceIds->InsertNextId(pts[faceVerts[4]]);
 567               faceIds->InsertNextId(pts[faceVerts[5]]);
 568               aCellType = VTK_POLYGON;
 569               numFacePts = 6;
 570             }
 571             input->GetCellNeighbors(cellId, faceIds, cellIds);
 572 #ifdef SHOW_COINCIDING_3D_PAL21924
 573             int nbNeighbors = 0;
 574             for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 575               if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 576             }
 577             bool process = nbNeighbors <= 0;
 578 #else
 579             bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 580 #endif
 581             if ( process || myShowInside ||
 582                  (!allVisible && !cellVis[cellIds->GetId(0)]) )
 583             {
 584               for ( i=0; i < numFacePts; i++)
 585                 aNewPts[i] = pts[faceVerts[i]];
 586               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 587               if(myStoreMapping)
 588                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 589               outputCD->CopyData(cd,cellId,newCellId);
 590             }
 591             }
 592           break;
 593         }
 594         case VTK_PYRAMID: {
 595 #ifdef SHOW_COINCIDING_3D_PAL21924
 596           faceIdsTmp->Reset();
 597           for (int ai=0; ai<npts; ai++)
 598             faceIdsTmp->InsertNextId(pts[ai]);
 599           input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 600 #endif
 601           for (faceId = 0; faceId < 5; faceId++)
 602             {
 603             faceIds->Reset();
 604             faceVerts = vtkPyramid::GetFaceArray(faceId);
 605             faceIds->InsertNextId(pts[faceVerts[0]]);
 606             faceIds->InsertNextId(pts[faceVerts[1]]);
 607             faceIds->InsertNextId(pts[faceVerts[2]]);
 608             aCellType = VTK_TRIANGLE;
 609             numFacePts = 3;
 610             if (faceVerts[3] >= 0)
 611               {
 612               faceIds->InsertNextId(pts[faceVerts[3]]);
 613               aCellType = VTK_QUAD;
 614               numFacePts = 4;
 615               }
 616             input->GetCellNeighbors(cellId, faceIds, cellIds);
 617 #ifdef SHOW_COINCIDING_3D_PAL21924
 618             int nbNeighbors = 0;
 619             for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 620               if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 621             }
 622             bool process = nbNeighbors <= 0;
 623 #else
 624             bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 625 #endif
 626             if ( process || myShowInside ||
 627                  (!allVisible && !cellVis[cellIds->GetId(0)]) )
 628               {
 629               for ( i=0; i < numFacePts; i++)
 630                 aNewPts[i] = pts[faceVerts[i]];
 631               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 632               if(myStoreMapping)
 633                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 634               outputCD->CopyData(cd,cellId,newCellId);
 635               }
 636             }
 637           break;
 638         }
 639
 640 #if VTK_XVERSION > 50700
 641         case VTK_POLYHEDRON:
 642           {
 643             //MESSAGE("VTK_POLYHEDRON geometry filter");
 644             vtkIdType nFaces = 0;
 645             vtkIdType* ptIds = 0;
 646             int idp = 0;
 647             input->GetFaceStream(cellId, nFaces, ptIds);
 648 #ifdef SHOW_COINCIDING_3D_PAL21924
 649             faceIdsTmp->Reset();
 650             for (int ai=0; ai<npts; ai++)
 651               faceIdsTmp->InsertNextId(pts[ai]);
 652             input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 653 #endif
 654             for (faceId = 0; faceId < nFaces; faceId++)
 655               {
 656                 faceIds->Reset();
 657                 numFacePts = ptIds[idp];
 658                 //MESSAGE("numFacePts="<< numFacePts);
 659                 int pt0 = ++idp;
 660                 for (i = 0; i < numFacePts; i++)
 661                   {
 662                     //MESSAGE("ptIds[" << idp + i << "]=" << ptIds[idp + i]);
 663                     faceIds->InsertNextId(ptIds[idp + i]);
 664                   }
 665                 idp += numFacePts;
 666                 switch (numFacePts)
 667                   {
 668                   case 3:
 669                     aCellType = VTK_TRIANGLE;
 670                     break;
 671                   case 4:
 672                     aCellType = VTK_QUAD;
 673                     break;
 674                   default:
 675                     aCellType = VTK_POLYGON;
 676                     break;
 677                   }
 678                 // TODO understand and fix display of several polyhedrons
 679                 input->GetCellNeighbors(cellId, faceIds, cellIds);
 680 #ifdef SHOW_COINCIDING_3D_PAL21924
 681                 int nbNeighbors = 0;
 682                 for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 683                   if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 684                 }
 685                 bool process = nbNeighbors <= 0;
 686 #else
 687                 bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 688 #endif
 689                 if (process || myShowInside
 690                     || (!allVisible && !cellVis[cellIds->GetId(0)]))
 691                   {
 692                     for (i = 0; i < numFacePts; i++)
 693                       aNewPts[i] = ptIds[pt0 + i];
 694                     newCellId = output->InsertNextCell(aCellType, numFacePts, aNewPts);
 695                     if (myStoreMapping)
 696                       InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 697                     outputCD->CopyData(cd, cellId, newCellId);
 698                   }
 699               }
 700             break;
 701           }
 702 #endif
 703         //Quadratic cells
 704         case VTK_QUADRATIC_EDGE:
 705         case VTK_QUADRATIC_TRIANGLE:
 706         case VTK_QUADRATIC_QUAD:
 707         case VTK_BIQUADRATIC_QUAD:
 708         case VTK_QUADRATIC_TETRA:
 709         case VTK_QUADRATIC_HEXAHEDRON:
 710         case VTK_TRIQUADRATIC_HEXAHEDRON:
 711         case VTK_QUADRATIC_WEDGE:
 712         case VTK_QUADRATIC_PYRAMID:
 713           if(!myIsWireframeMode) {
 714             input->GetCell(cellId,cell);
 715             vtkIdList *lpts = vtkIdList::New();
 716             vtkPoints *coords = vtkPoints::New();
 717             vtkIdList *cellIds = vtkIdList::New();
 718             vtkIdType newCellId;
 719
 720             if ( cell->GetCellDimension() == 1 ) {
 721               vtkIdType arcResult = -1;
 722               if(myIsBuildArc) {
 723                 arcResult = Build1DArc(cellId, input, output, pts, myMaxArcAngle);
 724                 newCellId = arcResult;
 725               }
 726
 727               if(!myIsBuildArc || arcResult == -1 ) {
 728                 aCellType = VTK_LINE;
 729                 numFacePts = 2;
 730                 cell->Triangulate(0,lpts,coords);
 731                 for (i=0; i < lpts->GetNumberOfIds(); i+=2) {
 732                   aNewPts[0] = lpts->GetId(i);
 733                   aNewPts[1] = lpts->GetId(i+1);
 734                   newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 735                   if(myStoreMapping)
 736                     InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 737                   outputCD->CopyData(cd,cellId,newCellId);
 738                 }
 739               }
 740               else {
 741                 if(myStoreMapping)
 742                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 743                 outputCD->CopyData(cd,cellId,newCellId);
 744               }
 745             }
 746             else if ( cell->GetCellDimension() == 2 ) {
 747               if(!myIsBuildArc) {
 748                 aCellType = VTK_TRIANGLE;
 749                 numFacePts = 3;
 750                 cell->Triangulate(0,lpts,coords);
 751                 for (i=0; i < lpts->GetNumberOfIds(); i+=3) {
 752                   aNewPts[0] = lpts->GetId(i);
 753                   aNewPts[1] = lpts->GetId(i+1);
 754                   aNewPts[2] = lpts->GetId(i+2);
 755                   newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 756                   if(myStoreMapping)
 757                     InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 758                   outputCD->CopyData(cd,cellId,newCellId);
 759                 }
 760               }
 761               else{
 762                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds,true);
 763               }
 764             }
 765             else //3D nonlinear cell
 766             {
 767 #ifdef SHOW_COINCIDING_3D_PAL21924
 768               faceIdsTmp->Reset();
 769               int npts1 = 0;
 770               switch (aCellType ){
 771               case VTK_QUADRATIC_TETRA:         npts1 = 4; break;
 772               case VTK_QUADRATIC_HEXAHEDRON:    npts1 = 8; break;
 773               case VTK_TRIQUADRATIC_HEXAHEDRON: npts1 = 8; break;
 774               case VTK_QUADRATIC_WEDGE:         npts1 = 6; break;
 775               case VTK_QUADRATIC_PYRAMID:       npts1 = 5; break;
 776               }
 777               if ( npts1 > 0 ) {
 778                 for (int ai=0; ai<npts; ai++)
 779                   faceIdsTmp->InsertNextId(pts[ai]);
 780                 input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 781               }
 782 #endif
 783               aCellType = VTK_TRIANGLE;
 784               numFacePts = 3;
 785               for (int j=0; j < cell->GetNumberOfFaces(); j++){
 786                 vtkCell *face = cell->GetFace(j);
 787                 input->GetCellNeighbors(cellId, face->PointIds, cellIds);
 788 #ifdef SHOW_COINCIDING_3D_PAL21924
 789                 int nbNeighbors = 0;
 790                 for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 791                   if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 792                 }
 793                 bool process = nbNeighbors <= 0;
 794 #else
 795                 bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 796 #endif
 797                 if ( process || myShowInside ) {
 798                   face->Triangulate(0,lpts,coords);
 799                   for (i=0; i < lpts->GetNumberOfIds(); i+=3) {
 800                     aNewPts[0] = lpts->GetId(i);
 801                     aNewPts[1] = lpts->GetId(i+1);
 802                     aNewPts[2] = lpts->GetId(i+2);
 803                     newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 804                     if(myStoreMapping)
 805                       InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 806                     outputCD->CopyData(cd,cellId,newCellId);
 807                   }
 808                 }
 809               }
 810             } //3d nonlinear cell
 811             cellIds->Delete();
 812             coords->Delete();
 813             lpts->Delete();
 814             break;
 815           } else { // wireframe
 816             switch(aCellType) {
 817             case VTK_QUADRATIC_EDGE: {
 818               vtkIdType arcResult =-1;
 819               if(myIsBuildArc) {
 820                arcResult = Build1DArc(cellId, input, output, pts,myMaxArcAngle);
 821                newCellId = arcResult;
 822               }
 823               if(!myIsBuildArc || arcResult == -1) {
 824                 aCellType = VTK_POLY_LINE;
 825                 numFacePts = 3;
 826
 827                 aNewPts[0] = pts[0];
 828                 aNewPts[2] = pts[1];
 829                 aNewPts[1] = pts[2];
 830
 831                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 832               }
 833
 834               if(myStoreMapping)
 835                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 836
 837               outputCD->CopyData(cd,cellId,newCellId);
 838               break;
 839             }
 840             case VTK_QUADRATIC_TRIANGLE: {
 841               if(!myIsBuildArc) {
 842                 aCellType = VTK_POLYGON;
 843                 numFacePts = 6;
 844
 845                 aNewPts[0] = pts[0];
 846                 aNewPts[1] = pts[3];
 847                 aNewPts[2] = pts[1];
 848                 aNewPts[3] = pts[4];
 849                 aNewPts[4] = pts[2];
 850                 aNewPts[5] = pts[5];
 851
 852                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 853                 if(myStoreMapping)
 854                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 855
 856                 outputCD->CopyData(cd,cellId,newCellId);
 857               }
 858               else
 859                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
 860               break;
 861             }
 862             case VTK_BIQUADRATIC_QUAD:
 863             case VTK_QUADRATIC_QUAD: {
 864               if(!myIsBuildArc) {
 865                 aCellType = VTK_POLYGON;
 866                 numFacePts = 8;
 867
 868                 aNewPts[0] = pts[0];
 869                 aNewPts[1] = pts[4];
 870                 aNewPts[2] = pts[1];
 871                 aNewPts[3] = pts[5];
 872                 aNewPts[4] = pts[2];
 873                 aNewPts[5] = pts[6];
 874                 aNewPts[6] = pts[3];
 875                 aNewPts[7] = pts[7];
 876
 877                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 878                 if(myStoreMapping)
 879                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 880
 881                 outputCD->CopyData(cd,cellId,newCellId);
 882               }
 883               else
 884                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
 885               break;
 886             }
 887             case VTK_QUADRATIC_TETRA: {
 888               aCellType = VTK_POLYGON;
 889               numFacePts = 6;
 890
 891               //---------------------------------------------------------------
 892               aNewPts[0] = pts[0];
 893               aNewPts[1] = pts[4];
 894               aNewPts[2] = pts[1];
 895               aNewPts[3] = pts[5];
 896               aNewPts[4] = pts[2];
 897               aNewPts[5] = pts[6];
 898
 899               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 900               if(myStoreMapping)
 901                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 902
 903               outputCD->CopyData(cd,cellId,newCellId);
 904
 905               //---------------------------------------------------------------
 906               aNewPts[0] = pts[0];
 907               aNewPts[1] = pts[7];
 908               aNewPts[2] = pts[3];
 909               aNewPts[3] = pts[8];
 910               aNewPts[4] = pts[1];
 911               aNewPts[5] = pts[4];
 912
 913               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 914               if(myStoreMapping)
 915                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 916
 917               outputCD->CopyData(cd,cellId,newCellId);
 918
 919               //---------------------------------------------------------------
 920               aNewPts[0] = pts[1];
 921               aNewPts[1] = pts[8];
 922               aNewPts[2] = pts[3];
 923               aNewPts[3] = pts[9];
 924               aNewPts[4] = pts[2];
 925               aNewPts[5] = pts[5];
 926
 927               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 928               if(myStoreMapping)
 929                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 930
 931               outputCD->CopyData(cd,cellId,newCellId);
 932
 933               //---------------------------------------------------------------
 934               aNewPts[0] = pts[2];
 935               aNewPts[1] = pts[9];
 936               aNewPts[2] = pts[3];
 937               aNewPts[3] = pts[7];
 938               aNewPts[4] = pts[0];
 939               aNewPts[5] = pts[6];
 940
 941               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 942               if(myStoreMapping)
 943                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 944
 945               outputCD->CopyData(cd,cellId,newCellId);
 946
 947               break;
 948             }
 949             case VTK_QUADRATIC_WEDGE: {
 950               aCellType = VTK_POLYGON;
 951               numFacePts = 6;
 952               //---------------------------------------------------------------
 953               //Face 1
 954               aNewPts[0] = pts[0];
 955               aNewPts[1] = pts[6];
 956               aNewPts[2] = pts[1];
 957               aNewPts[3] = pts[7];
 958               aNewPts[4] = pts[2];
 959               aNewPts[5] = pts[8];
 960               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 961               if(myStoreMapping)
 962                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 963               outputCD->CopyData(cd,cellId,newCellId);
 964
 965               //---------------------------------------------------------------
 966               //Face 2
 967               aNewPts[0] = pts[3];
 968               aNewPts[1] = pts[9];
 969               aNewPts[2] = pts[4];
 970               aNewPts[3] = pts[10];
 971               aNewPts[4] = pts[5];
 972               aNewPts[5] = pts[11];
 973               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 974               if(myStoreMapping)
 975                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 976               outputCD->CopyData(cd,cellId,newCellId);
 977
 978               //---------------------------------------------------------------
 979               //Face 3
 980               numFacePts = 8;
 981               aNewPts[0] = pts[0];
 982               aNewPts[1] = pts[8];
 983               aNewPts[2] = pts[2];
 984               aNewPts[3] = pts[14];
 985               aNewPts[4] = pts[5];
 986               aNewPts[5] = pts[11];
 987               aNewPts[6] = pts[3];
 988               aNewPts[7] = pts[12];
 989               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 990               if(myStoreMapping)
 991                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 992               outputCD->CopyData(cd,cellId,newCellId);
 993
 994               //---------------------------------------------------------------
 995               //Face 4
 996               aNewPts[0] = pts[1];
 997               aNewPts[1] = pts[13];
 998               aNewPts[2] = pts[4];
 999               aNewPts[3] = pts[10];
1000               aNewPts[4] = pts[5];
1001               aNewPts[5] = pts[14];
1002               aNewPts[6] = pts[2];
1003               aNewPts[7] = pts[7];
1004               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1005               if(myStoreMapping)
1006                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1007               outputCD->CopyData(cd,cellId,newCellId);
1008
1009               //---------------------------------------------------------------
1010               //Face 5
1011               aNewPts[0] = pts[0];
1012               aNewPts[1] = pts[12];
1013               aNewPts[2] = pts[3];
1014               aNewPts[3] = pts[9];
1015               aNewPts[4] = pts[4];
1016               aNewPts[5] = pts[13];
1017               aNewPts[6] = pts[1];
1018               aNewPts[7] = pts[6];
1019               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1020               if(myStoreMapping)
1021                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1022               outputCD->CopyData(cd,cellId,newCellId);
1023               break;
1024             }
1025             case VTK_TRIQUADRATIC_HEXAHEDRON:
1026             case VTK_QUADRATIC_HEXAHEDRON: {
1027               aCellType = VTK_POLYGON;
1028               numFacePts = 8;
1029
1030               //---------------------------------------------------------------
1031               aNewPts[0] = pts[0];
1032               aNewPts[1] = pts[8];
1033               aNewPts[2] = pts[1];
1034               aNewPts[3] = pts[17];
1035               aNewPts[4] = pts[5];
1036               aNewPts[5] = pts[12];
1037               aNewPts[6] = pts[4];
1038               aNewPts[7] = pts[16];
1039
1040               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1041               if(myStoreMapping)
1042                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1043
1044               outputCD->CopyData(cd,cellId,newCellId);
1045
1046               //---------------------------------------------------------------
1047               aNewPts[0] = pts[1];
1048               aNewPts[1] = pts[9];
1049               aNewPts[2] = pts[2];
1050               aNewPts[3] = pts[18];
1051               aNewPts[4] = pts[6];
1052               aNewPts[5] = pts[13];
1053               aNewPts[6] = pts[5];
1054               aNewPts[7] = pts[17];
1055
1056               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1057               if(myStoreMapping)
1058                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1059
1060               outputCD->CopyData(cd,cellId,newCellId);
1061
1062               //---------------------------------------------------------------
1063               aNewPts[0] = pts[2];
1064               aNewPts[1] = pts[10];
1065               aNewPts[2] = pts[3];
1066               aNewPts[3] = pts[19];
1067               aNewPts[4] = pts[7];
1068               aNewPts[5] = pts[14];
1069               aNewPts[6] = pts[6];
1070               aNewPts[7] = pts[18];
1071
1072               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1073               if(myStoreMapping)
1074                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1075
1076               outputCD->CopyData(cd,cellId,newCellId);
1077
1078               //---------------------------------------------------------------
1079               aNewPts[0] = pts[3];
1080               aNewPts[1] = pts[11];
1081               aNewPts[2] = pts[0];
1082               aNewPts[3] = pts[16];
1083               aNewPts[4] = pts[4];
1084               aNewPts[5] = pts[15];
1085               aNewPts[6] = pts[7];
1086               aNewPts[7] = pts[19];
1087
1088               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1089               if(myStoreMapping)
1090                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1091
1092               outputCD->CopyData(cd,cellId,newCellId);
1093
1094               //---------------------------------------------------------------
1095               aNewPts[0] = pts[0];
1096               aNewPts[1] = pts[8];
1097               aNewPts[2] = pts[1];
1098               aNewPts[3] = pts[9];
1099               aNewPts[4] = pts[2];
1100               aNewPts[5] = pts[10];
1101               aNewPts[6] = pts[3];
1102               aNewPts[7] = pts[11];
1103
1104               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1105               if(myStoreMapping)
1106                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1107
1108               outputCD->CopyData(cd,cellId,newCellId);
1109
1110               //---------------------------------------------------------------
1111               aNewPts[0] = pts[4];
1112               aNewPts[1] = pts[12];
1113               aNewPts[2] = pts[5];
1114               aNewPts[3] = pts[13];
1115               aNewPts[4] = pts[6];
1116               aNewPts[5] = pts[14];
1117               aNewPts[6] = pts[7];
1118               aNewPts[7] = pts[15];
1119
1120               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1121               if(myStoreMapping)
1122                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1123
1124               outputCD->CopyData(cd,cellId,newCellId);
1125
1126               break;
1127             }
1128             case VTK_QUADRATIC_PYRAMID: {
1129               aCellType = VTK_POLYGON;
1130               numFacePts = 6;
1131
1132               //---------------------------------------------------------------
1133               aNewPts[0] = pts[0];
1134               aNewPts[1] = pts[8];
1135               aNewPts[2] = pts[3];
1136               aNewPts[3] = pts[12];
1137               aNewPts[4] = pts[4];
1138               aNewPts[5] = pts[9];
1139
1140               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1141               if(myStoreMapping)
1142                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1143
1144               outputCD->CopyData(cd,cellId,newCellId);
1145
1146               //---------------------------------------------------------------
1147               aNewPts[0] = pts[0];
1148               aNewPts[1] = pts[9];
1149               aNewPts[2] = pts[4];
1150               aNewPts[3] = pts[10];
1151               aNewPts[4] = pts[1];
1152               aNewPts[5] = pts[5];
1153
1154               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1155               if(myStoreMapping)
1156                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1157
1158               outputCD->CopyData(cd,cellId,newCellId);
1159
1160               //---------------------------------------------------------------
1161               aNewPts[0] = pts[1];
1162               aNewPts[1] = pts[10];
1163               aNewPts[2] = pts[4];
1164               aNewPts[3] = pts[11];
1165               aNewPts[4] = pts[2];
1166               aNewPts[5] = pts[6];
1167
1168               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1169               if(myStoreMapping)
1170                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1171
1172               outputCD->CopyData(cd,cellId,newCellId);
1173
1174               //---------------------------------------------------------------
1175               aNewPts[0] = pts[2];
1176               aNewPts[1] = pts[11];
1177               aNewPts[2] = pts[4];
1178               aNewPts[3] = pts[12];
1179               aNewPts[4] = pts[3];
1180               aNewPts[5] = pts[7];
1181
1182               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1183               if(myStoreMapping)
1184                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1185
1186               outputCD->CopyData(cd,cellId,newCellId);
1187
1188               //---------------------------------------------------------------
1189               numFacePts = 8;
1190               aNewPts[0] = pts[0];
1191               aNewPts[1] = pts[5];
1192               aNewPts[2] = pts[1];
1193               aNewPts[3] = pts[6];
1194               aNewPts[4] = pts[2];
1195               aNewPts[5] = pts[7];
1196               aNewPts[6] = pts[3];
1197               aNewPts[7] = pts[8];
1198
1199               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1200               if(myStoreMapping)
1201                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1202
1203               outputCD->CopyData(cd,cellId,newCellId);
1204
1205               break;
1206             } // case VTK_QUADRATIC_PYRAMID:
1207             } // switch by type
1208           } // end WIREFRAME
1209           break;
1210         } //switch
1211       } //if visible
1212     } //for all cells
1213
1214   output->Squeeze();
1215
1216   vtkDebugMacro(<<"Extracted " << input->GetNumberOfPoints() << " points,"
1217                 << output->GetNumberOfCells() << " cells.");
1218
1219   cell->Delete();
1220
1221   cellIds->Delete();
1222   faceIds->Delete();
1223   cellIdsTmp->Delete();
1224   faceIdsTmp->Delete();
1225
1226   if ( cellVis )
1227   {
1228     delete [] cellVis;
1229   }
1230
1231   // to sort myVTK2ObjIds vector by cell dimension (ascending)
1232   if( myStoreMapping && !aDimension2VTK2ObjIds.empty() )
1233   {
1234     myVTK2ObjIds.clear();
1235     for( vtkIdType aDimension = 0; aDimension <= 2; aDimension++ )
1236     {
1237       TMapOfVectorId::const_iterator anIter = aDimension2VTK2ObjIds.find( aDimension );
1238       if( anIter != aDimension2VTK2ObjIds.end() )
1239       {
1240         const TVectorId& aCellIds = anIter->second;
1241         TVectorId::const_iterator anIdIter, anIdIterEnd = aCellIds.end();
1242         for( anIdIter = aCellIds.begin(); anIdIter != anIdIterEnd; anIdIter++ )
1243         {
1244           const vtkIdType aCellId = *anIdIter;
1245           myVTK2ObjIds.push_back( aCellId );
1246         }
1247       }
1248     }
1249   }
1250
1251   return 1;
1252 }
1253
1254 void
1255 VTKViewer_GeometryFilter
1256 ::InsertId( const vtkIdType theCellId,
1257             const vtkIdType theCellType,
1258             TVectorId& theVTK2ObjIds,
1259             TMapOfVectorId& theDimension2VTK2ObjIds )
1260 {
1261   theVTK2ObjIds.push_back( theCellId );
1262
1263   int aDimension = 0;
1264   switch( theCellType )
1265   {
1266     case VTK_VERTEX:
1267     case VTK_POLY_VERTEX:
1268       aDimension = 0;
1269       break;
1270     case VTK_LINE:
1271     case VTK_POLY_LINE:
1272       aDimension = 1;
1273       break;
1274     case VTK_TRIANGLE:
1275     case VTK_TRIANGLE_STRIP:
1276     case VTK_POLYGON:
1277     case VTK_PIXEL:
1278     case VTK_QUAD:
1279       aDimension = 2;
1280       break;
1281   }
1282
1283   TVectorId& aCellIds = theDimension2VTK2ObjIds[ aDimension ];
1284   aCellIds.push_back( theCellId );
1285 }
1286
1287 void
1288 VTKViewer_GeometryFilter
1289 ::SetInside(int theShowInside)
1290 {
1291   if(myShowInside == theShowInside)
1292     return;
1293
1294   myShowInside = theShowInside;
1295   this->Modified();
1296 }
1297
1298 int
1299 VTKViewer_GeometryFilter
1300 ::GetInside()
1301 {
1302   return myShowInside;
1303 }
1304
1305
1306 void
1307 VTKViewer_GeometryFilter
1308 ::SetWireframeMode(int theIsWireframeMode)
1309 {
1310   if(myIsWireframeMode == theIsWireframeMode)
1311     return;
1312
1313   myIsWireframeMode = theIsWireframeMode;
1314   this->Modified();
1315 }
1316
1317 int
1318 VTKViewer_GeometryFilter
1319 ::GetWireframeMode()
1320 {
1321   return myIsWireframeMode;
1322 }
1323
1324
1325 void
1326 VTKViewer_GeometryFilter
1327 ::SetStoreMapping(int theStoreMapping)
1328 {
1329   if(myStoreMapping == theStoreMapping)
1330     return;
1331
1332   myStoreMapping = theStoreMapping;
1333   this->Modified();
1334 }
1335
1336 int
1337 VTKViewer_GeometryFilter
1338 ::GetStoreMapping()
1339 {
1340   return myStoreMapping;
1341 }
1342
1343
1344 vtkIdType VTKViewer_GeometryFilter::GetElemObjId( int theVtkID )
1345 {
1346   if( theVtkID < 0 || theVtkID >= (int)myVTK2ObjIds.size() )
1347     return -1;
1348   return myVTK2ObjIds[theVtkID];
1349 }
1350
1351
1352 void VTKViewer_GeometryFilter::BuildArcedPolygon(vtkIdType cellId,
1353                                                  vtkUnstructuredGrid* input,
1354                                                  vtkPolyData *output,
1355                                                  TMapOfVectorId& theDimension2VTK2ObjIds,
1356                                                  bool triangulate){
1357   vtkIdType aCellType = VTK_POLYGON;
1358   vtkIdType *aNewPoints = NULL;
1359   vtkIdType aNbPoints = 0;
1360   vtkIdType newCellId;
1361
1362   //Input and output cell data
1363   vtkCellData *cd = input->GetCellData();
1364   vtkCellData *outputCD = output->GetCellData();
1365
1366   //Input and output scalars on point data
1367   vtkDataArray* inputScalars = input->GetPointData()->GetScalars();
1368   vtkDataArray* outputScalars = output->GetPointData()->GetScalars();
1369
1370   std::vector<vtkPoints*> aCollection;
1371   std::vector< std::vector<double> > aScalarCollection;
1372
1373   vtkCell* aCell = input->GetCell(cellId);
1374   switch(aCell->GetCellType()) {
1375     case VTK_QUADRATIC_TRIANGLE:
1376     {
1377       //Get All points from input cell
1378       Pnt P0 = CreatePnt( aCell, inputScalars, 0 );
1379       Pnt P1 = CreatePnt( aCell, inputScalars, 1 );
1380       Pnt P2 = CreatePnt( aCell, inputScalars, 2 );
1381       Pnt P3 = CreatePnt( aCell, inputScalars, 3 );
1382       Pnt P4 = CreatePnt( aCell, inputScalars, 4 );
1383       Pnt P5 = CreatePnt( aCell, inputScalars, 5 );
1384
1385       VTKViewer_ArcBuilder aBuilder1(P0,P3,P1,myMaxArcAngle); //Build arc using 0, 3 and 1 points
1386 #ifdef __MYDEBUG__
1387       cout << "Quadrangle arc 1 " << ( aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1388 #endif
1389
1390       VTKViewer_ArcBuilder aBuilder2(P1,P4,P2,myMaxArcAngle); //Build arc using 1, 4 and 2 points
1391 #ifdef __MYDEBUG__
1392       cout << "Quadrangle arc 2 " << ( aBuilder2.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1393 #endif
1394
1395       VTKViewer_ArcBuilder aBuilder3(P2,P5,P0,myMaxArcAngle); //Build arc using 2, 5 and 0 points
1396 #ifdef __MYDEBUG__
1397       cout << "Quadrangle arc 3 " << ( aBuilder3.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1398 #endif
1399
1400       aCollection.push_back(aBuilder1.GetPoints());
1401       aCollection.push_back(aBuilder2.GetPoints());
1402       aCollection.push_back(aBuilder3.GetPoints());
1403
1404       aScalarCollection.push_back(aBuilder1.GetScalarValues());
1405       aScalarCollection.push_back(aBuilder2.GetScalarValues());
1406       aScalarCollection.push_back(aBuilder3.GetScalarValues());
1407       break;
1408     }
1409     case VTK_BIQUADRATIC_QUAD:
1410     case VTK_QUADRATIC_QUAD:
1411     {
1412       //Get All points from input cell
1413       Pnt P0 = CreatePnt( aCell, inputScalars, 0 );
1414       Pnt P1 = CreatePnt( aCell, inputScalars, 1 );
1415       Pnt P2 = CreatePnt( aCell, inputScalars, 2 );
1416       Pnt P3 = CreatePnt( aCell, inputScalars, 3 );
1417       Pnt P4 = CreatePnt( aCell, inputScalars, 4 );
1418       Pnt P5 = CreatePnt( aCell, inputScalars, 5 );
1419       Pnt P6 = CreatePnt( aCell, inputScalars, 6 );
1420       Pnt P7 = CreatePnt( aCell, inputScalars, 7 );
1421
1422       VTKViewer_ArcBuilder aBuilder1(P0,P4,P1,myMaxArcAngle); //Build arc using 0, 4 and 1 points
1423 #ifdef __MYDEBUG__
1424       cout << "Quadrangle arc 1 " << ( aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1425 #endif
1426
1427       VTKViewer_ArcBuilder aBuilder2(P1,P5,P2,myMaxArcAngle); //Build arc using 1, 5 and 2 points
1428 #ifdef __MYDEBUG__
1429       cout << "Quadrangle arc 2 " << ( aBuilder2.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1430 #endif
1431
1432       VTKViewer_ArcBuilder aBuilder3(P2,P6,P3,myMaxArcAngle); //Build arc using 2, 6 and 3 points
1433 #ifdef __MYDEBUG__
1434       cout << "Quadrangle arc 3 " << ( aBuilder3.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1435 #endif
1436
1437       VTKViewer_ArcBuilder aBuilder4(P3,P7,P0,myMaxArcAngle); //Build arc using 3, 7 and 0 points
1438 #ifdef __MYDEBUG__
1439       cout << "Quadrangle arc 4 " << ( aBuilder4.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1440 #endif
1441
1442       aCollection.push_back(aBuilder1.GetPoints());
1443       aCollection.push_back(aBuilder2.GetPoints());
1444       aCollection.push_back(aBuilder3.GetPoints());
1445       aCollection.push_back(aBuilder4.GetPoints());
1446
1447       aScalarCollection.push_back(aBuilder1.GetScalarValues());
1448       aScalarCollection.push_back(aBuilder2.GetScalarValues());
1449       aScalarCollection.push_back(aBuilder3.GetScalarValues());
1450       aScalarCollection.push_back(aBuilder4.GetScalarValues());
1451       break;
1452     }
1453     default: //Unsupported cell type
1454       return;
1455   }
1456
1457   if(triangulate){
1458     const vtkIdType numFacePts = 3;
1459     vtkIdList *pts = vtkIdList::New();
1460     vtkPoints *coords = vtkPoints::New();
1461     aCellType = VTK_TRIANGLE;
1462     vtkIdType aNewPts[numFacePts];
1463     vtkIdType aTriangleId;
1464
1465     vtkPolygon *aPlg = vtkPolygon::New();
1466     std::map<int, double> aPntId2ScalarValue;
1467     aNbPoints = MergevtkPoints(aCollection, aScalarCollection, aPlg->GetPoints(), aPntId2ScalarValue, aNewPoints);
1468     aPlg->GetPointIds()->SetNumberOfIds(aNbPoints);
1469
1470     for(vtkIdType i = 0; i < aNbPoints;i++) {
1471       aPlg->GetPointIds()->SetId(i, aNewPoints[i]);
1472     }
1473
1474     aPlg->Triangulate(0,pts,coords);
1475
1476     for (vtkIdType i=0; i < pts->GetNumberOfIds(); i+=3) {
1477       aNewPts[0] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i));
1478       aNewPts[1] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+1));
1479       aNewPts[2] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+2));
1480
1481       if(outputScalars) {
1482         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i)]);
1483         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i+1)]);
1484         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i+2)]);
1485       }
1486
1487       aTriangleId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1488
1489       if(myStoreMapping)
1490         InsertId( cellId, aCellType, myVTK2ObjIds, theDimension2VTK2ObjIds );
1491       outputCD->CopyData(cd,cellId,aTriangleId);
1492     }
1493     pts->Delete();
1494     coords->Delete();
1495     aPlg->Delete();
1496   }
1497   else {
1498     std::map<int, double> aPntId2ScalarValue;
1499     aNbPoints = MergevtkPoints(aCollection, aScalarCollection, output->GetPoints(), aPntId2ScalarValue, aNewPoints);
1500     if(outputScalars)
1501       for(vtkIdType i = 0; i < aNbPoints; i++)
1502         outputScalars->InsertNextTuple1(aPntId2ScalarValue[aNewPoints[i]]);
1503     newCellId = output->InsertNextCell(aCellType,aNbPoints,aNewPoints);
1504     outputCD->CopyData(cd,cellId,newCellId);
1505
1506     if(myStoreMapping)
1507       InsertId( cellId, aCellType, myVTK2ObjIds, theDimension2VTK2ObjIds );
1508   }
1509
1510   if (aNewPoints)
1511     delete [] aNewPoints;
1512 }
1513
1514
1515 void VTKViewer_GeometryFilter::SetQuadraticArcMode(bool theFlag)
1516 {
1517   if(myIsBuildArc != theFlag) {
1518     myIsBuildArc = theFlag;
1519     this->Modified();
1520   }
1521 }
1522 bool VTKViewer_GeometryFilter::GetQuadraticArcMode() const
1523 {
1524   return myIsBuildArc;
1525 }
1526
1527 void VTKViewer_GeometryFilter::SetQuadraticArcAngle(double theMaxAngle)
1528 {
1529   if(myMaxArcAngle != theMaxAngle) {
1530     myMaxArcAngle = theMaxAngle;
1531     this->Modified();
1532   }
1533 }
1534
1535 double VTKViewer_GeometryFilter:: GetQuadraticArcAngle() const
1536 {
1537   return myMaxArcAngle;
1538 }
1539
1540 int VTKViewer_GeometryFilter::GetAppendCoincident3D() const {
1541 // VSR 26/10/2012: see description of SHOW_COINCIDING_3D_PAL20314
1542 // in the top of this file
1543 #ifdef SHOW_COINCIDING_3D_PAL20314
1544   return myAppendCoincident3D;
1545 #else
1546   return false;
1547 #endif
1548 }
1549
1550 void VTKViewer_GeometryFilter::SetAppendCoincident3D(int theFlag) {
1551   if(myAppendCoincident3D != theFlag){
1552     myAppendCoincident3D = theFlag;
1553     this->Modified();
1554   }
1555 }