src/VTKViewer/VTKViewer_GeometryFilter.cxx

   1 // Copyright (C) 2007-2015  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, or (at your option) any later version.
  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 #include <vtkVersion.h>
  52
  53 #include <algorithm>
  54 #include <iterator>
  55 #include <vector>
  56 #include <map>
  57 #include <set>
  58
  59 #include "utilities.h"
  60
  61 #if defined __GNUC__
  62   #if __GNUC__ == 2
  63     #define __GNUC_2__
  64   #endif
  65 #endif
  66
  67 #define VTK_XVERSION (VTK_MAJOR_VERSION*10000+VTK_MINOR_VERSION*100+VTK_BUILD_VERSION)
  68
  69 //#define __MYDEBUG__
  70 //#define USE_ROBUST_TRIANGULATION
  71
  72 ///////////////////////////////////////////////////////////////////////////////////////////////
  73 // VSR 26/10/2012: fix of regression (issue 21924) - increased memory consumption
  74 // for displaying of 3d elements, introduced by fix for issue 20314.
  75 // ...
  76 // The macro SHOW_COINCIDING_3D_PAL20314, when defined, allows correct visualization of
  77 // coincident 3d elements but causes substantial increasing of memory consumption, as all 3d
  78 // elements are always shown, even if they are totally covered by surrounding faces.
  79 // If this macro is not defined (commented), the behaviour is defined by another macro -
  80 // SHOW_COINCIDING_3D_PAL21924, as follows:
  81 // - If SHOW_COINCIDING_3D_PAL21924 is defined, an alternative solution for computing
  82 //   visibility of 3d elements is used; this solution allows to fix problem with visibility
  83 //   of coinciding 3d elements in most cases (though some cases might not work), while not
  84 //   causing significant increasing of memory consumption.
  85 // - If SHOW_COINCIDING_3D_PAL21924 is not defined (commented), coinciding 3d elements are
  86 //   not shown at all (this corresponds to the state before issue 20314 fixing).
  87 ///////////////////////////////////////////////////////////////////////////////////////////////
  88 //#define SHOW_COINCIDING_3D_PAL20314
  89 #ifndef SHOW_COINCIDING_3D_PAL20314
  90 #define SHOW_COINCIDING_3D_PAL21924
  91 #endif
  92 ///////////////////////////////////////////////////////////////////////////////////////////////
  93
  94 vtkStandardNewMacro(VTKViewer_GeometryFilter);
  95
  96 VTKViewer_GeometryFilter
  97 ::VTKViewer_GeometryFilter():
  98   myShowInside(0),
  99   myStoreMapping(0),
 100   myAppendCoincident3D(0),
 101   myIsWireframeMode(0),
 102   myIsBuildArc(false),
 103   myMaxArcAngle(2)
 104 {}
 105
 106
 107 VTKViewer_GeometryFilter
 108 ::~VTKViewer_GeometryFilter()
 109 {}
 110
 111
 112 int
 113 VTKViewer_GeometryFilter
 114 ::RequestData(
 115   vtkInformation *request,
 116   vtkInformationVector **inputVector,
 117   vtkInformationVector *outputVector)
 118 {
 119   // get the info objects
 120   vtkInformation *inInfo = inputVector[0]->GetInformationObject(0);
 121   vtkInformation *outInfo = outputVector->GetInformationObject(0);
 122
 123   // get the input and ouptut
 124   vtkDataSet *input = vtkDataSet::SafeDownCast(
 125     inInfo->Get(vtkDataObject::DATA_OBJECT()));
 126   vtkPolyData *output = vtkPolyData::SafeDownCast(
 127     outInfo->Get(vtkDataObject::DATA_OBJECT()));
 128
 129   vtkIdType numCells=input->GetNumberOfCells();
 130
 131   if (numCells == 0)
 132     {
 133       return 0;
 134     }
 135
 136   if (input->GetDataObjectType() == VTK_UNSTRUCTURED_GRID){
 137     return this->UnstructuredGridExecute(input, output, outInfo);
 138   }else
 139     return Superclass::RequestData(request,inputVector,outputVector);
 140
 141   return 1;
 142 }
 143
 144
 145 int
 146 VTKViewer_GeometryFilter
 147 ::UnstructuredGridExecute(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 boundaries)
 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_BIQUADRATIC_TRIANGLE:
 707         case VTK_QUADRATIC_QUAD:
 708         case VTK_BIQUADRATIC_QUAD:
 709         case VTK_QUADRATIC_POLYGON:
 710         case VTK_QUADRATIC_TETRA:
 711         case VTK_QUADRATIC_HEXAHEDRON:
 712         case VTK_TRIQUADRATIC_HEXAHEDRON:
 713         case VTK_QUADRATIC_WEDGE:
 714         case VTK_QUADRATIC_PYRAMID:
 715           if(!myIsWireframeMode) {
 716             input->GetCell(cellId,cell);
 717             vtkIdList *lpts = vtkIdList::New();
 718             vtkPoints *coords = vtkPoints::New();
 719             vtkIdList *cellIds = vtkIdList::New();
 720             vtkIdType newCellId;
 721
 722             if ( cell->GetCellDimension() == 1 ) {
 723               vtkIdType arcResult = -1;
 724               if(myIsBuildArc) {
 725                 arcResult = Build1DArc(cellId, input, output, pts, myMaxArcAngle);
 726                 newCellId = arcResult;
 727               }
 728
 729               if(!myIsBuildArc || arcResult == -1 ) {
 730                 aCellType = VTK_LINE;
 731                 numFacePts = 2;
 732                 cell->Triangulate(0,lpts,coords);
 733                 for (i=0; i < lpts->GetNumberOfIds(); i+=2) {
 734                   aNewPts[0] = lpts->GetId(i);
 735                   aNewPts[1] = lpts->GetId(i+1);
 736                   newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 737                   if(myStoreMapping)
 738                     InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 739                   outputCD->CopyData(cd,cellId,newCellId);
 740                 }
 741               }
 742               else {
 743                 if(myStoreMapping)
 744                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 745                 outputCD->CopyData(cd,cellId,newCellId);
 746               }
 747             }
 748             else if ( cell->GetCellDimension() == 2 ) {
 749               if(!myIsBuildArc) {
 750                 aCellType = VTK_TRIANGLE;
 751                 numFacePts = 3;
 752                 cell->Triangulate(0,lpts,coords);
 753                 for (i=0; i < lpts->GetNumberOfIds(); i+=3) {
 754                   aNewPts[0] = lpts->GetId(i);
 755                   aNewPts[1] = lpts->GetId(i+1);
 756                   aNewPts[2] = lpts->GetId(i+2);
 757                   newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 758                   if(myStoreMapping)
 759                     InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 760                   outputCD->CopyData(cd,cellId,newCellId);
 761                 }
 762               }
 763               else{
 764                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds,true);
 765               }
 766             }
 767             else //3D nonlinear cell
 768             {
 769 #ifdef SHOW_COINCIDING_3D_PAL21924
 770               faceIdsTmp->Reset();
 771               int npts1 = 0;
 772               switch (aCellType ){
 773               case VTK_QUADRATIC_TETRA:         npts1 = 4; break;
 774               case VTK_QUADRATIC_HEXAHEDRON:    npts1 = 8; break;
 775               case VTK_TRIQUADRATIC_HEXAHEDRON: npts1 = 8; break;
 776               case VTK_QUADRATIC_WEDGE:         npts1 = 6; break;
 777               case VTK_QUADRATIC_PYRAMID:       npts1 = 5; break;
 778               }
 779               if ( npts1 > 0 ) {
 780                 for (int ai=0; ai<npts; ai++)
 781                   faceIdsTmp->InsertNextId(pts[ai]);
 782                 input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 783               }
 784 #endif
 785               aCellType = VTK_TRIANGLE;
 786               numFacePts = 3;
 787               for (int j=0; j < cell->GetNumberOfFaces(); j++){
 788                 vtkCell *face = cell->GetFace(j);
 789                 input->GetCellNeighbors(cellId, face->PointIds, cellIds);
 790 #ifdef SHOW_COINCIDING_3D_PAL21924
 791                 int nbNeighbors = 0;
 792                 for(int ai=0;ai<cellIds->GetNumberOfIds();ai++) {
 793                   if (cellIdsTmp->IsId(cellIds->GetId(ai)) == -1) nbNeighbors++;
 794                 }
 795                 bool process = nbNeighbors <= 0;
 796 #else
 797                 bool process = cellIds->GetNumberOfIds() <= 0 || GetAppendCoincident3D();
 798 #endif
 799                 if ( process || myShowInside ) {
 800                   face->Triangulate(0,lpts,coords);
 801                   for (i=0; i < lpts->GetNumberOfIds(); i+=3) {
 802                     aNewPts[0] = lpts->GetId(i);
 803                     aNewPts[1] = lpts->GetId(i+1);
 804                     aNewPts[2] = lpts->GetId(i+2);
 805                     newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 806                     if(myStoreMapping)
 807                       InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 808                     outputCD->CopyData(cd,cellId,newCellId);
 809                   }
 810                 }
 811               }
 812             } //3d nonlinear cell
 813             cellIds->Delete();
 814             coords->Delete();
 815             lpts->Delete();
 816             break;
 817           } else { // wireframe
 818             switch(aCellType) {
 819             case VTK_QUADRATIC_EDGE: {
 820               vtkIdType arcResult =-1;
 821               if(myIsBuildArc) {
 822                arcResult = Build1DArc(cellId, input, output, pts,myMaxArcAngle);
 823                newCellId = arcResult;
 824               }
 825               if(!myIsBuildArc || arcResult == -1) {
 826                 aCellType = VTK_POLY_LINE;
 827                 numFacePts = 3;
 828
 829                 aNewPts[0] = pts[0];
 830                 aNewPts[2] = pts[1];
 831                 aNewPts[1] = pts[2];
 832
 833                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 834               }
 835
 836               if(myStoreMapping)
 837                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 838
 839               outputCD->CopyData(cd,cellId,newCellId);
 840               break;
 841             }
 842             case VTK_QUADRATIC_TRIANGLE:
 843             case VTK_BIQUADRATIC_TRIANGLE: {
 844               if(!myIsBuildArc) {
 845                 aCellType = VTK_POLYGON;
 846                 numFacePts = 6;
 847
 848                 aNewPts[0] = pts[0];
 849                 aNewPts[1] = pts[3];
 850                 aNewPts[2] = pts[1];
 851                 aNewPts[3] = pts[4];
 852                 aNewPts[4] = pts[2];
 853                 aNewPts[5] = pts[5];
 854
 855                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 856                 if(myStoreMapping)
 857                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 858
 859                 outputCD->CopyData(cd,cellId,newCellId);
 860               }
 861               else
 862                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
 863               break;
 864             }
 865             case VTK_QUADRATIC_QUAD:
 866             case VTK_BIQUADRATIC_QUAD: {
 867               if(!myIsBuildArc) {
 868                 aCellType = VTK_POLYGON;
 869                 numFacePts = 8;
 870
 871                 aNewPts[0] = pts[0];
 872                 aNewPts[1] = pts[4];
 873                 aNewPts[2] = pts[1];
 874                 aNewPts[3] = pts[5];
 875                 aNewPts[4] = pts[2];
 876                 aNewPts[5] = pts[6];
 877                 aNewPts[6] = pts[3];
 878                 aNewPts[7] = pts[7];
 879
 880                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 881                 if(myStoreMapping)
 882                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 883
 884                 outputCD->CopyData(cd,cellId,newCellId);
 885               }
 886               else
 887                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
 888               break;
 889             }
 890             case VTK_QUADRATIC_POLYGON: {
 891               if(!myIsBuildArc)
 892               {
 893                 aCellType = VTK_POLYGON;
 894
 895                 for ( i = 0; i < npts/2; ++i )
 896                 {
 897                   aNewPts[i*2  ] = pts[i];
 898                   aNewPts[i*2+1] = pts[i+npts/2];
 899                 }
 900                 newCellId = output->InsertNextCell(aCellType,npts,aNewPts);
 901                 if(myStoreMapping)
 902                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 903
 904                 outputCD->CopyData(cd,cellId,newCellId);
 905               }
 906               else
 907                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
 908               break;
 909             }
 910             case VTK_QUADRATIC_TETRA: {
 911               aCellType = VTK_POLYGON;
 912               numFacePts = 6;
 913
 914               //---------------------------------------------------------------
 915               aNewPts[0] = pts[0];
 916               aNewPts[1] = pts[4];
 917               aNewPts[2] = pts[1];
 918               aNewPts[3] = pts[5];
 919               aNewPts[4] = pts[2];
 920               aNewPts[5] = pts[6];
 921
 922               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 923               if(myStoreMapping)
 924                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 925
 926               outputCD->CopyData(cd,cellId,newCellId);
 927
 928               //---------------------------------------------------------------
 929               aNewPts[0] = pts[0];
 930               aNewPts[1] = pts[7];
 931               aNewPts[2] = pts[3];
 932               aNewPts[3] = pts[8];
 933               aNewPts[4] = pts[1];
 934               aNewPts[5] = pts[4];
 935
 936               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 937               if(myStoreMapping)
 938                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 939
 940               outputCD->CopyData(cd,cellId,newCellId);
 941
 942               //---------------------------------------------------------------
 943               aNewPts[0] = pts[1];
 944               aNewPts[1] = pts[8];
 945               aNewPts[2] = pts[3];
 946               aNewPts[3] = pts[9];
 947               aNewPts[4] = pts[2];
 948               aNewPts[5] = pts[5];
 949
 950               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 951               if(myStoreMapping)
 952                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 953
 954               outputCD->CopyData(cd,cellId,newCellId);
 955
 956               //---------------------------------------------------------------
 957               aNewPts[0] = pts[2];
 958               aNewPts[1] = pts[9];
 959               aNewPts[2] = pts[3];
 960               aNewPts[3] = pts[7];
 961               aNewPts[4] = pts[0];
 962               aNewPts[5] = pts[6];
 963
 964               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 965               if(myStoreMapping)
 966                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 967
 968               outputCD->CopyData(cd,cellId,newCellId);
 969
 970               break;
 971             }
 972             case VTK_QUADRATIC_WEDGE: {
 973               aCellType = VTK_POLYGON;
 974               numFacePts = 6;
 975               //---------------------------------------------------------------
 976               //Face 1
 977               aNewPts[0] = pts[0];
 978               aNewPts[1] = pts[6];
 979               aNewPts[2] = pts[1];
 980               aNewPts[3] = pts[7];
 981               aNewPts[4] = pts[2];
 982               aNewPts[5] = pts[8];
 983               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 984               if(myStoreMapping)
 985                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 986               outputCD->CopyData(cd,cellId,newCellId);
 987
 988               //---------------------------------------------------------------
 989               //Face 2
 990               aNewPts[0] = pts[3];
 991               aNewPts[1] = pts[9];
 992               aNewPts[2] = pts[4];
 993               aNewPts[3] = pts[10];
 994               aNewPts[4] = pts[5];
 995               aNewPts[5] = pts[11];
 996               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 997               if(myStoreMapping)
 998                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 999               outputCD->CopyData(cd,cellId,newCellId);
1000
1001               //---------------------------------------------------------------
1002               //Face 3
1003               numFacePts = 8;
1004               aNewPts[0] = pts[0];
1005               aNewPts[1] = pts[8];
1006               aNewPts[2] = pts[2];
1007               aNewPts[3] = pts[14];
1008               aNewPts[4] = pts[5];
1009               aNewPts[5] = pts[11];
1010               aNewPts[6] = pts[3];
1011               aNewPts[7] = pts[12];
1012               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1013               if(myStoreMapping)
1014                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1015               outputCD->CopyData(cd,cellId,newCellId);
1016
1017               //---------------------------------------------------------------
1018               //Face 4
1019               aNewPts[0] = pts[1];
1020               aNewPts[1] = pts[13];
1021               aNewPts[2] = pts[4];
1022               aNewPts[3] = pts[10];
1023               aNewPts[4] = pts[5];
1024               aNewPts[5] = pts[14];
1025               aNewPts[6] = pts[2];
1026               aNewPts[7] = pts[7];
1027               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1028               if(myStoreMapping)
1029                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1030               outputCD->CopyData(cd,cellId,newCellId);
1031
1032               //---------------------------------------------------------------
1033               //Face 5
1034               aNewPts[0] = pts[0];
1035               aNewPts[1] = pts[12];
1036               aNewPts[2] = pts[3];
1037               aNewPts[3] = pts[9];
1038               aNewPts[4] = pts[4];
1039               aNewPts[5] = pts[13];
1040               aNewPts[6] = pts[1];
1041               aNewPts[7] = pts[6];
1042               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1043               if(myStoreMapping)
1044                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1045               outputCD->CopyData(cd,cellId,newCellId);
1046               break;
1047             }
1048             case VTK_TRIQUADRATIC_HEXAHEDRON:
1049             case VTK_QUADRATIC_HEXAHEDRON: {
1050               aCellType = VTK_POLYGON;
1051               numFacePts = 8;
1052
1053               //---------------------------------------------------------------
1054               aNewPts[0] = pts[0];
1055               aNewPts[1] = pts[8];
1056               aNewPts[2] = pts[1];
1057               aNewPts[3] = pts[17];
1058               aNewPts[4] = pts[5];
1059               aNewPts[5] = pts[12];
1060               aNewPts[6] = pts[4];
1061               aNewPts[7] = pts[16];
1062
1063               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1064               if(myStoreMapping)
1065                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1066
1067               outputCD->CopyData(cd,cellId,newCellId);
1068
1069               //---------------------------------------------------------------
1070               aNewPts[0] = pts[1];
1071               aNewPts[1] = pts[9];
1072               aNewPts[2] = pts[2];
1073               aNewPts[3] = pts[18];
1074               aNewPts[4] = pts[6];
1075               aNewPts[5] = pts[13];
1076               aNewPts[6] = pts[5];
1077               aNewPts[7] = pts[17];
1078
1079               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1080               if(myStoreMapping)
1081                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1082
1083               outputCD->CopyData(cd,cellId,newCellId);
1084
1085               //---------------------------------------------------------------
1086               aNewPts[0] = pts[2];
1087               aNewPts[1] = pts[10];
1088               aNewPts[2] = pts[3];
1089               aNewPts[3] = pts[19];
1090               aNewPts[4] = pts[7];
1091               aNewPts[5] = pts[14];
1092               aNewPts[6] = pts[6];
1093               aNewPts[7] = pts[18];
1094
1095               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1096               if(myStoreMapping)
1097                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1098
1099               outputCD->CopyData(cd,cellId,newCellId);
1100
1101               //---------------------------------------------------------------
1102               aNewPts[0] = pts[3];
1103               aNewPts[1] = pts[11];
1104               aNewPts[2] = pts[0];
1105               aNewPts[3] = pts[16];
1106               aNewPts[4] = pts[4];
1107               aNewPts[5] = pts[15];
1108               aNewPts[6] = pts[7];
1109               aNewPts[7] = pts[19];
1110
1111               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1112               if(myStoreMapping)
1113                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1114
1115               outputCD->CopyData(cd,cellId,newCellId);
1116
1117               //---------------------------------------------------------------
1118               aNewPts[0] = pts[0];
1119               aNewPts[1] = pts[8];
1120               aNewPts[2] = pts[1];
1121               aNewPts[3] = pts[9];
1122               aNewPts[4] = pts[2];
1123               aNewPts[5] = pts[10];
1124               aNewPts[6] = pts[3];
1125               aNewPts[7] = pts[11];
1126
1127               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1128               if(myStoreMapping)
1129                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1130
1131               outputCD->CopyData(cd,cellId,newCellId);
1132
1133               //---------------------------------------------------------------
1134               aNewPts[0] = pts[4];
1135               aNewPts[1] = pts[12];
1136               aNewPts[2] = pts[5];
1137               aNewPts[3] = pts[13];
1138               aNewPts[4] = pts[6];
1139               aNewPts[5] = pts[14];
1140               aNewPts[6] = pts[7];
1141               aNewPts[7] = pts[15];
1142
1143               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1144               if(myStoreMapping)
1145                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1146
1147               outputCD->CopyData(cd,cellId,newCellId);
1148
1149               break;
1150             }
1151             case VTK_QUADRATIC_PYRAMID: {
1152               aCellType = VTK_POLYGON;
1153               numFacePts = 6;
1154
1155               //---------------------------------------------------------------
1156               aNewPts[0] = pts[0];
1157               aNewPts[1] = pts[8];
1158               aNewPts[2] = pts[3];
1159               aNewPts[3] = pts[12];
1160               aNewPts[4] = pts[4];
1161               aNewPts[5] = pts[9];
1162
1163               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1164               if(myStoreMapping)
1165                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1166
1167               outputCD->CopyData(cd,cellId,newCellId);
1168
1169               //---------------------------------------------------------------
1170               aNewPts[0] = pts[0];
1171               aNewPts[1] = pts[9];
1172               aNewPts[2] = pts[4];
1173               aNewPts[3] = pts[10];
1174               aNewPts[4] = pts[1];
1175               aNewPts[5] = pts[5];
1176
1177               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1178               if(myStoreMapping)
1179                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1180
1181               outputCD->CopyData(cd,cellId,newCellId);
1182
1183               //---------------------------------------------------------------
1184               aNewPts[0] = pts[1];
1185               aNewPts[1] = pts[10];
1186               aNewPts[2] = pts[4];
1187               aNewPts[3] = pts[11];
1188               aNewPts[4] = pts[2];
1189               aNewPts[5] = pts[6];
1190
1191               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1192               if(myStoreMapping)
1193                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1194
1195               outputCD->CopyData(cd,cellId,newCellId);
1196
1197               //---------------------------------------------------------------
1198               aNewPts[0] = pts[2];
1199               aNewPts[1] = pts[11];
1200               aNewPts[2] = pts[4];
1201               aNewPts[3] = pts[12];
1202               aNewPts[4] = pts[3];
1203               aNewPts[5] = pts[7];
1204
1205               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1206               if(myStoreMapping)
1207                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1208
1209               outputCD->CopyData(cd,cellId,newCellId);
1210
1211               //---------------------------------------------------------------
1212               numFacePts = 8;
1213               aNewPts[0] = pts[0];
1214               aNewPts[1] = pts[5];
1215               aNewPts[2] = pts[1];
1216               aNewPts[3] = pts[6];
1217               aNewPts[4] = pts[2];
1218               aNewPts[5] = pts[7];
1219               aNewPts[6] = pts[3];
1220               aNewPts[7] = pts[8];
1221
1222               newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1223               if(myStoreMapping)
1224                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1225
1226               outputCD->CopyData(cd,cellId,newCellId);
1227
1228               break;
1229             } // case VTK_QUADRATIC_PYRAMID:
1230             } // switch by type
1231           } // end WIREFRAME
1232           break;
1233         } //switch
1234       } //if visible
1235     } //for all cells
1236
1237   output->Squeeze();
1238
1239   vtkDebugMacro(<<"Extracted " << input->GetNumberOfPoints() << " points,"
1240                 << output->GetNumberOfCells() << " cells.");
1241
1242   cell->Delete();
1243
1244   cellIds->Delete();
1245   faceIds->Delete();
1246   cellIdsTmp->Delete();
1247   faceIdsTmp->Delete();
1248
1249   if ( cellVis )
1250   {
1251     delete [] cellVis;
1252   }
1253
1254   // to sort myVTK2ObjIds vector by cell dimension (ascending)
1255   if( myStoreMapping && !aDimension2VTK2ObjIds.empty() )
1256   {
1257     myVTK2ObjIds.clear();
1258     for( vtkIdType aDimension = 0; aDimension <= 2; aDimension++ )
1259     {
1260       TMapOfVectorId::const_iterator anIter = aDimension2VTK2ObjIds.find( aDimension );
1261       if( anIter != aDimension2VTK2ObjIds.end() )
1262       {
1263         const TVectorId& aCellIds = anIter->second;
1264         TVectorId::const_iterator anIdIter, anIdIterEnd = aCellIds.end();
1265         for( anIdIter = aCellIds.begin(); anIdIter != anIdIterEnd; anIdIter++ )
1266         {
1267           const vtkIdType aCellId = *anIdIter;
1268           myVTK2ObjIds.push_back( aCellId );
1269         }
1270       }
1271     }
1272   }
1273
1274   return 1;
1275 }
1276
1277 void
1278 VTKViewer_GeometryFilter
1279 ::InsertId( const vtkIdType theCellId,
1280             const vtkIdType theCellType,
1281             TVectorId& theVTK2ObjIds,
1282             TMapOfVectorId& theDimension2VTK2ObjIds )
1283 {
1284   theVTK2ObjIds.push_back( theCellId );
1285
1286   int aDimension = 0;
1287   switch( theCellType )
1288   {
1289     case VTK_VERTEX:
1290     case VTK_POLY_VERTEX:
1291       aDimension = 0;
1292       break;
1293     case VTK_LINE:
1294     case VTK_POLY_LINE:
1295       aDimension = 1;
1296       break;
1297     case VTK_TRIANGLE:
1298     case VTK_TRIANGLE_STRIP:
1299     case VTK_POLYGON:
1300     case VTK_PIXEL:
1301     case VTK_QUAD:
1302       aDimension = 2;
1303       break;
1304   }
1305
1306   TVectorId& aCellIds = theDimension2VTK2ObjIds[ aDimension ];
1307   aCellIds.push_back( theCellId );
1308 }
1309
1310 void
1311 VTKViewer_GeometryFilter
1312 ::SetInside(int theShowInside)
1313 {
1314   if(myShowInside == theShowInside)
1315     return;
1316
1317   myShowInside = theShowInside;
1318   this->Modified();
1319 }
1320
1321 int
1322 VTKViewer_GeometryFilter
1323 ::GetInside()
1324 {
1325   return myShowInside;
1326 }
1327
1328
1329 void
1330 VTKViewer_GeometryFilter
1331 ::SetWireframeMode(int theIsWireframeMode)
1332 {
1333   if(myIsWireframeMode == theIsWireframeMode)
1334     return;
1335
1336   myIsWireframeMode = theIsWireframeMode;
1337   this->Modified();
1338 }
1339
1340 int
1341 VTKViewer_GeometryFilter
1342 ::GetWireframeMode()
1343 {
1344   return myIsWireframeMode;
1345 }
1346
1347
1348 void
1349 VTKViewer_GeometryFilter
1350 ::SetStoreMapping(int theStoreMapping)
1351 {
1352   if(myStoreMapping == theStoreMapping)
1353     return;
1354
1355   myStoreMapping = theStoreMapping;
1356   this->Modified();
1357 }
1358
1359 int
1360 VTKViewer_GeometryFilter
1361 ::GetStoreMapping()
1362 {
1363   return myStoreMapping;
1364 }
1365
1366
1367 vtkIdType VTKViewer_GeometryFilter::GetElemObjId( int theVtkID )
1368 {
1369   if( theVtkID < 0 || theVtkID >= (int)myVTK2ObjIds.size() )
1370     return -1;
1371   return myVTK2ObjIds[theVtkID];
1372 }
1373
1374
1375 void VTKViewer_GeometryFilter::BuildArcedPolygon(vtkIdType cellId,
1376                                                  vtkUnstructuredGrid* input,
1377                                                  vtkPolyData *output,
1378                                                  TMapOfVectorId& theDimension2VTK2ObjIds,
1379                                                  bool triangulate)
1380 {
1381   vtkIdType aCellType = VTK_POLYGON;
1382   vtkIdType *aNewPoints = NULL;
1383   vtkIdType aNbPoints = 0;
1384   vtkIdType newCellId;
1385
1386   //Input and output cell data
1387   vtkCellData *cd = input->GetCellData();
1388   vtkCellData *outputCD = output->GetCellData();
1389
1390   //Input and output scalars on point data
1391   vtkDataArray* inputScalars = input->GetPointData()->GetScalars();
1392   vtkDataArray* outputScalars = output->GetPointData()->GetScalars();
1393
1394   std::vector< vtkSmartPointer<vtkPoints> > aCollection;
1395   std::vector< std::vector<double> > aScalarCollection;
1396
1397   vtkCell* aCell = input->GetCell(cellId);
1398   switch(aCell->GetCellType()) {
1399     case VTK_QUADRATIC_TRIANGLE:
1400     case VTK_BIQUADRATIC_TRIANGLE:
1401     {
1402       //Get All points from input cell
1403       Pnt P0 = CreatePnt( aCell, inputScalars, 0 );
1404       Pnt P1 = CreatePnt( aCell, inputScalars, 1 );
1405       Pnt P2 = CreatePnt( aCell, inputScalars, 2 );
1406       Pnt P3 = CreatePnt( aCell, inputScalars, 3 );
1407       Pnt P4 = CreatePnt( aCell, inputScalars, 4 );
1408       Pnt P5 = CreatePnt( aCell, inputScalars, 5 );
1409
1410       VTKViewer_ArcBuilder aBuilder1(P0,P3,P1,myMaxArcAngle); //Build arc using 0, 3 and 1 points
1411 #ifdef __MYDEBUG__
1412       cout << "Quadrangle arc 1 " << ( aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1413 #endif
1414
1415       VTKViewer_ArcBuilder aBuilder2(P1,P4,P2,myMaxArcAngle); //Build arc using 1, 4 and 2 points
1416 #ifdef __MYDEBUG__
1417       cout << "Quadrangle arc 2 " << ( aBuilder2.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1418 #endif
1419
1420       VTKViewer_ArcBuilder aBuilder3(P2,P5,P0,myMaxArcAngle); //Build arc using 2, 5 and 0 points
1421 #ifdef __MYDEBUG__
1422       cout << "Quadrangle arc 3 " << ( aBuilder3.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1423 #endif
1424
1425       aCollection.push_back(aBuilder1.GetPoints());
1426       aCollection.push_back(aBuilder2.GetPoints());
1427       aCollection.push_back(aBuilder3.GetPoints());
1428
1429       aScalarCollection.push_back(aBuilder1.GetScalarValues());
1430       aScalarCollection.push_back(aBuilder2.GetScalarValues());
1431       aScalarCollection.push_back(aBuilder3.GetScalarValues());
1432       break;
1433     }
1434     case VTK_QUADRATIC_QUAD:
1435     case VTK_BIQUADRATIC_QUAD:
1436     {
1437       //Get All points from input cell
1438       Pnt P0 = CreatePnt( aCell, inputScalars, 0 );
1439       Pnt P1 = CreatePnt( aCell, inputScalars, 1 );
1440       Pnt P2 = CreatePnt( aCell, inputScalars, 2 );
1441       Pnt P3 = CreatePnt( aCell, inputScalars, 3 );
1442       Pnt P4 = CreatePnt( aCell, inputScalars, 4 );
1443       Pnt P5 = CreatePnt( aCell, inputScalars, 5 );
1444       Pnt P6 = CreatePnt( aCell, inputScalars, 6 );
1445       Pnt P7 = CreatePnt( aCell, inputScalars, 7 );
1446
1447       VTKViewer_ArcBuilder aBuilder1(P0,P4,P1,myMaxArcAngle); //Build arc using 0, 4 and 1 points
1448 #ifdef __MYDEBUG__
1449       cout << "Quadrangle arc 1 " << ( aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1450 #endif
1451
1452       VTKViewer_ArcBuilder aBuilder2(P1,P5,P2,myMaxArcAngle); //Build arc using 1, 5 and 2 points
1453 #ifdef __MYDEBUG__
1454       cout << "Quadrangle arc 2 " << ( aBuilder2.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1455 #endif
1456
1457       VTKViewer_ArcBuilder aBuilder3(P2,P6,P3,myMaxArcAngle); //Build arc using 2, 6 and 3 points
1458 #ifdef __MYDEBUG__
1459       cout << "Quadrangle arc 3 " << ( aBuilder3.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1460 #endif
1461
1462       VTKViewer_ArcBuilder aBuilder4(P3,P7,P0,myMaxArcAngle); //Build arc using 3, 7 and 0 points
1463 #ifdef __MYDEBUG__
1464       cout << "Quadrangle arc 4 " << ( aBuilder4.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1465 #endif
1466
1467       aCollection.push_back(aBuilder1.GetPoints());
1468       aCollection.push_back(aBuilder2.GetPoints());
1469       aCollection.push_back(aBuilder3.GetPoints());
1470       aCollection.push_back(aBuilder4.GetPoints());
1471
1472       aScalarCollection.push_back(aBuilder1.GetScalarValues());
1473       aScalarCollection.push_back(aBuilder2.GetScalarValues());
1474       aScalarCollection.push_back(aBuilder3.GetScalarValues());
1475       aScalarCollection.push_back(aBuilder4.GetScalarValues());
1476       break;
1477     }
1478     case VTK_QUADRATIC_POLYGON:
1479     {
1480       int nbP = aCell->GetNumberOfPoints();
1481       std::vector< Pnt > pVec( nbP + 2 );
1482
1483       for ( int i = 0; i < nbP/2; ++i )
1484       {
1485         pVec[i*2 + 0] = CreatePnt( aCell, inputScalars, i );
1486         pVec[i*2 + 1] = CreatePnt( aCell, inputScalars, i + nbP/2 );
1487       }
1488       pVec[ nbP   ] = pVec[ 0 ];
1489       pVec[ nbP+1 ] = pVec[ 1 ];
1490
1491       for ( int i = 0; i < nbP; i += 2 )
1492       {
1493         VTKViewer_ArcBuilder aBuilder( pVec[i], pVec[i+1], pVec[i+2], myMaxArcAngle );
1494         aCollection.push_back( aBuilder.GetPoints() );
1495         aScalarCollection.push_back( aBuilder.GetScalarValues() );
1496       }
1497       break;
1498     }
1499     default: //Unsupported cell type
1500       return;
1501   }
1502
1503   if(triangulate){
1504     const vtkIdType numFacePts = 3;
1505     vtkIdList *pts = vtkIdList::New();
1506     vtkPoints *coords = vtkPoints::New();
1507     aCellType = VTK_TRIANGLE;
1508     vtkIdType aNewPts[numFacePts];
1509     vtkIdType aTriangleId;
1510
1511     vtkPolygon *aPlg = vtkPolygon::New();
1512     std::map<int, double> aPntId2ScalarValue;
1513     aNbPoints = MergevtkPoints(aCollection, aScalarCollection, aPlg->GetPoints(), aPntId2ScalarValue, aNewPoints);
1514     aPlg->GetPointIds()->SetNumberOfIds(aNbPoints);
1515
1516     for(vtkIdType i = 0; i < aNbPoints;i++) {
1517       aPlg->GetPointIds()->SetId(i, aNewPoints[i]);
1518     }
1519
1520     aPlg->Triangulate(0,pts,coords);
1521
1522     for (vtkIdType i=0; i < pts->GetNumberOfIds(); i+=3) {
1523       aNewPts[0] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i));
1524       aNewPts[1] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+1));
1525       aNewPts[2] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+2));
1526
1527       if(outputScalars) {
1528         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i)]);
1529         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i+1)]);
1530         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i+2)]);
1531       }
1532
1533       aTriangleId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1534
1535       if(myStoreMapping)
1536         InsertId( cellId, aCellType, myVTK2ObjIds, theDimension2VTK2ObjIds );
1537       outputCD->CopyData(cd,cellId,aTriangleId);
1538     }
1539     pts->Delete();
1540     coords->Delete();
1541     aPlg->Delete();
1542   }
1543   else {
1544     std::map<int, double> aPntId2ScalarValue;
1545     aNbPoints = MergevtkPoints(aCollection, aScalarCollection, output->GetPoints(), aPntId2ScalarValue, aNewPoints);
1546     if(outputScalars)
1547       for(vtkIdType i = 0; i < aNbPoints; i++)
1548         outputScalars->InsertNextTuple1(aPntId2ScalarValue[aNewPoints[i]]);
1549     newCellId = output->InsertNextCell(aCellType,aNbPoints,aNewPoints);
1550     outputCD->CopyData(cd,cellId,newCellId);
1551
1552     if(myStoreMapping)
1553       InsertId( cellId, aCellType, myVTK2ObjIds, theDimension2VTK2ObjIds );
1554   }
1555
1556   if (aNewPoints)
1557     delete [] aNewPoints;
1558 }
1559
1560
1561 void VTKViewer_GeometryFilter::SetQuadraticArcMode(bool theFlag)
1562 {
1563   if(myIsBuildArc != theFlag) {
1564     myIsBuildArc = theFlag;
1565     this->Modified();
1566   }
1567 }
1568 bool VTKViewer_GeometryFilter::GetQuadraticArcMode() const
1569 {
1570   return myIsBuildArc;
1571 }
1572
1573 void VTKViewer_GeometryFilter::SetQuadraticArcAngle(double theMaxAngle)
1574 {
1575   if(myMaxArcAngle != theMaxAngle) {
1576     myMaxArcAngle = theMaxAngle;
1577     this->Modified();
1578   }
1579 }
1580
1581 double VTKViewer_GeometryFilter:: GetQuadraticArcAngle() const
1582 {
1583   return myMaxArcAngle;
1584 }
1585
1586 int VTKViewer_GeometryFilter::GetAppendCoincident3D() const {
1587 // VSR 26/10/2012: see description of SHOW_COINCIDING_3D_PAL20314
1588 // in the top of this file
1589 #ifdef SHOW_COINCIDING_3D_PAL20314
1590   return myAppendCoincident3D;
1591 #else
1592   return false;
1593 #endif
1594 }
1595
1596 void VTKViewer_GeometryFilter::SetAppendCoincident3D(int theFlag) {
1597   if(myAppendCoincident3D != theFlag){
1598     myAppendCoincident3D = theFlag;
1599     this->Modified();
1600   }
1601 }