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