src/VTKViewer/VTKViewer_GeometryFilter.cxx

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