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_BIQUADRATIC_QUADRATIC_WEDGE:
 881         case VTK_QUADRATIC_PYRAMID:
 882
 883           if(!myIsWireframeMode)
 884           {
 885             input->GetCell(cellId,cell);
 886             vtkIdList *lpts = vtkIdList::New();
 887             vtkPoints *coords = vtkPoints::New();
 888             vtkIdList *cellIds = vtkIdList::New();
 889             vtkIdType newCellId;
 890
 891             if ( cell->GetCellDimension() == 1 ) {
 892               vtkIdType arcResult = -1;
 893               if(myIsBuildArc) {
 894                 arcResult = Build1DArc(cellId, input, output, pts, myMaxArcAngle);
 895                 newCellId = arcResult;
 896               }
 897
 898               if(!myIsBuildArc || arcResult == -1 ) {
 899                 aCellType = VTK_LINE;
 900                 numFacePts = 2;
 901                 cell->Triangulate(0,lpts,coords);
 902                 for (i=0; i < lpts->GetNumberOfIds(); i+=2) {
 903                   aNewPts[0] = lpts->GetId(i);
 904                   aNewPts[1] = lpts->GetId(i+1);
 905                   newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 906                   if(myStoreMapping)
 907                     InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 908                   outputCD->CopyData(cd,cellId,newCellId);
 909                 }
 910               }
 911               else {
 912                 if(myStoreMapping)
 913                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 914                 outputCD->CopyData(cd,cellId,newCellId);
 915               }
 916             }
 917             else if ( cell->GetCellDimension() == 2 ) {
 918               if(!myIsBuildArc) {
 919                 aCellType = VTK_TRIANGLE;
 920                 numFacePts = 3;
 921                 cell->Triangulate(0,lpts,coords);
 922                 for (i=0; i < lpts->GetNumberOfIds(); i+=3) {
 923                   aNewPts[0] = lpts->GetId(i);
 924                   aNewPts[1] = lpts->GetId(i+1);
 925                   aNewPts[2] = lpts->GetId(i+2);
 926                   newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 927                   if(myStoreMapping)
 928                     InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 929                   outputCD->CopyData(cd,cellId,newCellId);
 930                 }
 931               }
 932               else{
 933                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds,true);
 934               }
 935             }
 936             else //3D nonlinear cell
 937             {
 938 #ifdef SHOW_COINCIDING_3D_PAL21924
 939               if ( !myShowInside )
 940               {
 941                 int npts1 = 0;
 942                 switch (aCellType ){
 943                 case VTK_QUADRATIC_TETRA:             npts1 = 4; break;
 944                 case VTK_QUADRATIC_HEXAHEDRON:        npts1 = 8; break;
 945                 case VTK_TRIQUADRATIC_HEXAHEDRON:     npts1 = 8; break;
 946                 case VTK_QUADRATIC_WEDGE:             npts1 = 6; break;
 947                 case VTK_BIQUADRATIC_QUADRATIC_WEDGE: npts1 = 6; break;
 948                 case VTK_QUADRATIC_PYRAMID:           npts1 = 5; break;
 949                 }
 950                 faceIdsTmp->SetNumberOfIds( npts1 );
 951                 if ( npts1 > 0 ) {
 952                   for (int ai=0; ai<npts1; ai++)
 953                     faceIdsTmp->SetId( ai, pts[ai] );
 954                   input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
 955                 }
 956               }
 957 #endif
 958               aCellType = VTK_TRIANGLE;
 959               numFacePts = 3;
 960               int nbNeighbors = 0;
 961               for (int j=0; j < cell->GetNumberOfFaces(); j++)
 962               {
 963                 vtkCell *face = cell->GetFace(j);
 964                 if ( !myShowInside ) {
 965                   input->GetCellNeighbors(cellId, face->PointIds, cellIds);
 966                   nbNeighbors = cellIds->GetNumberOfIds() - cellIdsTmp->GetNumberOfIds();
 967                 }
 968 #ifdef SHOW_COINCIDING_3D_PAL21924
 969                 bool process = nbNeighbors <= 0;
 970 #else
 971                 bool process = nbNeighbors <= 0 || GetAppendCoincident3D();
 972 #endif
 973                 if ( process || myShowInside ) {
 974                   face->Triangulate(0,lpts,coords);
 975                   for (i=0; i < lpts->GetNumberOfIds(); i+=3) {
 976                     aNewPts[0] = lpts->GetId(i);
 977                     aNewPts[1] = lpts->GetId(i+1);
 978                     aNewPts[2] = lpts->GetId(i+2);
 979                     newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
 980                     if(myStoreMapping)
 981                       InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
 982                     outputCD->CopyData(cd,cellId,newCellId);
 983                   }
 984                 }
 985               }
 986             } //3d nonlinear cell
 987             cellIds->Delete();
 988             coords->Delete();
 989             lpts->Delete();
 990             break;
 991           }
 992           else { // wireframe
 993             switch(aCellType) {
 994             case VTK_QUADRATIC_EDGE:
 995             {
 996               vtkIdType arcResult =-1;
 997               if(myIsBuildArc) {
 998                arcResult = Build1DArc(cellId, input, output, pts,myMaxArcAngle);
 999                newCellId = arcResult;
1000               }
1001               if(!myIsBuildArc || arcResult == -1) {
1002                 aCellType = VTK_POLY_LINE;
1003                 numFacePts = 3;
1004
1005                 aNewPts[0] = pts[0];
1006                 aNewPts[2] = pts[1];
1007                 aNewPts[1] = pts[2];
1008
1009                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1010               }
1011
1012               if(myStoreMapping)
1013                 InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1014
1015               outputCD->CopyData(cd,cellId,newCellId);
1016               break;
1017             }
1018             case VTK_QUADRATIC_TRIANGLE:
1019             case VTK_BIQUADRATIC_TRIANGLE:
1020             {
1021               if(!myIsBuildArc) {
1022                 aCellType = VTK_POLYGON;
1023                 numFacePts = 6;
1024
1025                 aNewPts[0] = pts[0];
1026                 aNewPts[1] = pts[3];
1027                 aNewPts[2] = pts[1];
1028                 aNewPts[3] = pts[4];
1029                 aNewPts[4] = pts[2];
1030                 aNewPts[5] = pts[5];
1031
1032                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1033                 if(myStoreMapping)
1034                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1035
1036                 outputCD->CopyData(cd,cellId,newCellId);
1037               }
1038               else
1039                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
1040               break;
1041             }
1042             case VTK_QUADRATIC_QUAD:
1043             case VTK_BIQUADRATIC_QUAD:
1044             {
1045               if(!myIsBuildArc) {
1046                 aCellType = VTK_POLYGON;
1047                 numFacePts = 8;
1048
1049                 aNewPts[0] = pts[0];
1050                 aNewPts[1] = pts[4];
1051                 aNewPts[2] = pts[1];
1052                 aNewPts[3] = pts[5];
1053                 aNewPts[4] = pts[2];
1054                 aNewPts[5] = pts[6];
1055                 aNewPts[6] = pts[3];
1056                 aNewPts[7] = pts[7];
1057
1058                 newCellId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1059                 if(myStoreMapping)
1060                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1061
1062                 outputCD->CopyData(cd,cellId,newCellId);
1063               }
1064               else
1065                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
1066               break;
1067             }
1068             case VTK_QUADRATIC_POLYGON:
1069             {
1070               if(!myIsBuildArc)
1071               {
1072                 aCellType = VTK_POLYGON;
1073
1074                 for ( i = 0; i < npts/2; ++i )
1075                 {
1076                   aNewPts[i*2  ] = pts[i];
1077                   aNewPts[i*2+1] = pts[i+npts/2];
1078                 }
1079                 newCellId = output->InsertNextCell(aCellType,npts,aNewPts);
1080                 if(myStoreMapping)
1081                   InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1082
1083                 outputCD->CopyData(cd,cellId,newCellId);
1084               }
1085               else
1086                 BuildArcedPolygon(cellId,input,output,aDimension2VTK2ObjIds);
1087               break;
1088             }
1089             case VTK_QUADRATIC_TETRA:
1090             case VTK_QUADRATIC_WEDGE:
1091             case VTK_BIQUADRATIC_QUADRATIC_WEDGE:
1092             case VTK_TRIQUADRATIC_HEXAHEDRON:
1093             case VTK_QUADRATIC_HEXAHEDRON:
1094             case VTK_QUADRATIC_PYRAMID:
1095             {
1096               aCellType = VTK_POLY_LINE;
1097               input->GetCell(cellId,cell);
1098               if ( myShowInside )
1099               {
1100                 int nbEdges = cell->GetNumberOfEdges();
1101                 for ( int edgeId = 0; edgeId < nbEdges; ++edgeId )
1102                 {
1103                   vtkCell * edge = cell->GetEdge( edgeId );
1104                   if ( toShowEdge( edge->GetPointId(0), edge->GetPointId(2), cellId, input ))
1105                   {
1106                     aNewPts[0] = edge->GetPointId(0);
1107                     aNewPts[1] = edge->GetPointId(2);
1108                     aNewPts[2] = edge->GetPointId(1);
1109                     newCellId = output->InsertNextCell( aCellType, 3, aNewPts );
1110                     if(myStoreMapping)
1111                       InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1112                     outputCD->CopyData(cd,cellId,newCellId);
1113                   }
1114                 }
1115               }
1116               else
1117               {
1118                 int nbCoincident = 0;
1119 #ifdef SHOW_COINCIDING_3D_PAL21924
1120                 int nbPnt = npts - cell->GetNumberOfEdges();
1121                 faceIdsTmp->SetNumberOfIds( nbPnt );
1122                 for ( int ai = 0; ai < nbPnt; ai++ )
1123                   faceIdsTmp->SetId( ai, pts[ai] );
1124                 input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
1125                 nbCoincident = cellIdsTmp->GetNumberOfIds();
1126 #endif
1127                 midPoints.clear();
1128                 int nbFaces = cell->GetNumberOfFaces();
1129                 for ( faceId = 0; faceId < nbFaces; faceId++ )
1130                 {
1131                   vtkCell * face = cell->GetFace( faceId );
1132                   input->GetCellNeighbors( cellId, face->GetPointIds(), cellIds );
1133                   int nbNeighbors = cellIds->GetNumberOfIds() - nbCoincident;
1134                   if ( nbNeighbors <= 0 )
1135                   {
1136                     int nbEdges = face->GetNumberOfPoints() / 2;
1137                     for ( int edgeId = 0; edgeId < nbEdges; ++edgeId )
1138                     {
1139                       vtkIdType p1 = ( edgeId );               // corner
1140                       vtkIdType p2 = ( edgeId + nbEdges );     // medium
1141                       vtkIdType p3 = ( edgeId + 1 ) % nbEdges; // next corner
1142                       faceIdsTmp->SetNumberOfIds( 2 );
1143                       faceIdsTmp->SetId( 0, face->GetPointId(p2) );
1144                       faceIdsTmp->SetId( 1, face->GetPointId(p1) );
1145                       input->GetCellNeighbors(cellId, faceIdsTmp, cellIdsTmp);
1146                       bool process;
1147                       switch ( cellIdsTmp->GetNumberOfIds() ) {
1148                       case 0: // the edge belong to this cell only
1149                         // avoid adding it when treating another face
1150                         process = midPoints.insert( face->GetPointId(p2) ).second; break;
1151                       case 1: // the edge is shared by two cells
1152                         process = ( cellIdsTmp->GetId(0) < cellId ); break;
1153                       default: // the edge is shared by >2 cells
1154                         process = ( cellIdsTmp->GetId(0) < cellId ); break;
1155                       }
1156                       if ( process )
1157                       {
1158                         aNewPts[0] = face->GetPointId( p1 );
1159                         aNewPts[1] = face->GetPointId( p2 );
1160                         aNewPts[2] = face->GetPointId( p3 );
1161                         newCellId = output->InsertNextCell( aCellType, 3, aNewPts );
1162                         if(myStoreMapping)
1163                           InsertId( cellId, aCellType, myVTK2ObjIds, aDimension2VTK2ObjIds );
1164                         outputCD->CopyData(cd,cellId,newCellId);
1165                       }
1166                     }
1167                   }
1168                 }
1169               }
1170               break;
1171             } // case of volumes in wireframe
1172             } // switch by quadratic type
1173           } // end WIREFRAME
1174           break;
1175         } //switch by type
1176
1177     } //if visible
1178   } //for all cells
1179
1180   output->Squeeze();
1181
1182   vtkDebugMacro(<<"Extracted " << input->GetNumberOfPoints() << " points,"
1183                 << output->GetNumberOfCells() << " cells.");
1184
1185   cell->Delete();
1186
1187   cellIds->Delete();
1188   faceIds->Delete();
1189   cellIdsTmp->Delete();
1190   faceIdsTmp->Delete();
1191
1192   if ( cellVis )
1193   {
1194     delete [] cellVis;
1195   }
1196
1197   if ( input->GetCellLinks() )
1198   {
1199     input->GetCellLinks()->Initialize(); // free memory
1200   }
1201
1202   // fill myVTK2ObjIds vector in ascending cell dimension order
1203   myVTK2ObjIds.clear();
1204   if( myStoreMapping && !aDimension2VTK2ObjIds.empty() )
1205   {
1206     size_t nbCells = ( aDimension2VTK2ObjIds[0].size() +
1207                        aDimension2VTK2ObjIds[1].size() +
1208                        aDimension2VTK2ObjIds[2].size() );
1209     if ( myVTK2ObjIds.capacity() > nbCells )
1210       TVectorId().swap( myVTK2ObjIds );
1211     myVTK2ObjIds.reserve( nbCells );
1212
1213     for( int aDimension = 0; aDimension <= 2; aDimension++ )
1214       if ( !aDimension2VTK2ObjIds[ aDimension ].empty() )
1215       {
1216         myVTK2ObjIds.insert( myVTK2ObjIds.end(),
1217                              aDimension2VTK2ObjIds[ aDimension ].begin(),
1218                              aDimension2VTK2ObjIds[ aDimension ].end() );
1219         TVectorId().swap( aDimension2VTK2ObjIds[ aDimension ]);
1220       }
1221   }
1222
1223   return 1;
1224 }
1225
1226 void
1227 VTKViewer_GeometryFilter
1228 ::InsertId( const vtkIdType theCellId,
1229             const vtkIdType theCellType,
1230             TVectorId& theVTK2ObjIds,
1231             TMapOfVectorId& theDimension2VTK2ObjIds )
1232 {
1233   //theVTK2ObjIds.push_back( theCellId );
1234
1235   int aDimension = 0;
1236   switch( theCellType )
1237   {
1238     case VTK_VERTEX:
1239     case VTK_POLY_VERTEX:
1240       aDimension = 0;
1241       break;
1242     case VTK_LINE:
1243     case VTK_POLY_LINE:
1244       aDimension = 1;
1245       break;
1246     case VTK_TRIANGLE:
1247     case VTK_TRIANGLE_STRIP:
1248     case VTK_POLYGON:
1249     case VTK_PIXEL:
1250     case VTK_QUAD:
1251       aDimension = 2;
1252       break;
1253   }
1254
1255   TVectorId& aCellIds = theDimension2VTK2ObjIds[ aDimension ];
1256   aCellIds.push_back( theCellId );
1257 }
1258
1259 void
1260 VTKViewer_GeometryFilter
1261 ::SetInside(int theShowInside)
1262 {
1263   if(myShowInside == theShowInside)
1264     return;
1265
1266   myShowInside = theShowInside;
1267   this->Modified();
1268 }
1269
1270 int
1271 VTKViewer_GeometryFilter
1272 ::GetInside()
1273 {
1274   return myShowInside;
1275 }
1276
1277
1278 void
1279 VTKViewer_GeometryFilter
1280 ::SetWireframeMode(int theIsWireframeMode)
1281 {
1282   if(myIsWireframeMode == theIsWireframeMode)
1283     return;
1284
1285   myIsWireframeMode = theIsWireframeMode;
1286   this->Modified();
1287 }
1288
1289 int
1290 VTKViewer_GeometryFilter
1291 ::GetWireframeMode()
1292 {
1293   return myIsWireframeMode;
1294 }
1295
1296
1297 void
1298 VTKViewer_GeometryFilter
1299 ::SetStoreMapping(int theStoreMapping)
1300 {
1301   if(myStoreMapping == theStoreMapping)
1302     return;
1303
1304   myStoreMapping = theStoreMapping;
1305   this->Modified();
1306 }
1307
1308 int
1309 VTKViewer_GeometryFilter
1310 ::GetStoreMapping()
1311 {
1312   return myStoreMapping;
1313 }
1314
1315
1316 vtkIdType VTKViewer_GeometryFilter::GetElemObjId( int theVtkID )
1317 {
1318   if( theVtkID < 0 || theVtkID >= (int)myVTK2ObjIds.size() )
1319     return -1;
1320   return myVTK2ObjIds[theVtkID];
1321 }
1322
1323
1324 void VTKViewer_GeometryFilter::BuildArcedPolygon(vtkIdType cellId,
1325                                                  vtkUnstructuredGrid* input,
1326                                                  vtkPolyData *output,
1327                                                  TMapOfVectorId& theDimension2VTK2ObjIds,
1328                                                  bool triangulate)
1329 {
1330   vtkIdType aCellType = VTK_POLYGON;
1331   vtkIdType *aNewPoints = NULL;
1332   vtkIdType aNbPoints = 0;
1333   vtkIdType newCellId;
1334
1335   //Input and output cell data
1336   vtkCellData *cd = input->GetCellData();
1337   vtkCellData *outputCD = output->GetCellData();
1338
1339   //Input and output scalars on point data
1340   vtkDataArray* inputScalars = input->GetPointData()->GetScalars();
1341   vtkDataArray* outputScalars = output->GetPointData()->GetScalars();
1342
1343   std::vector< vtkSmartPointer<vtkPoints> > aCollection;
1344   std::vector< std::vector<double> > aScalarCollection;
1345
1346   vtkCell* aCell = input->GetCell(cellId);
1347   switch(aCell->GetCellType()) {
1348     case VTK_QUADRATIC_TRIANGLE:
1349     case VTK_BIQUADRATIC_TRIANGLE:
1350     {
1351       //Get All points from input cell
1352       Pnt P0 = CreatePnt( aCell, inputScalars, 0 );
1353       Pnt P1 = CreatePnt( aCell, inputScalars, 1 );
1354       Pnt P2 = CreatePnt( aCell, inputScalars, 2 );
1355       Pnt P3 = CreatePnt( aCell, inputScalars, 3 );
1356       Pnt P4 = CreatePnt( aCell, inputScalars, 4 );
1357       Pnt P5 = CreatePnt( aCell, inputScalars, 5 );
1358
1359       VTKViewer_ArcBuilder aBuilder1(P0,P3,P1,myMaxArcAngle); //Build arc using 0, 3 and 1 points
1360 #ifdef __MYDEBUG__
1361       cout << "Quadrangle arc 1 " << ( aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1362 #endif
1363
1364       VTKViewer_ArcBuilder aBuilder2(P1,P4,P2,myMaxArcAngle); //Build arc using 1, 4 and 2 points
1365 #ifdef __MYDEBUG__
1366       cout << "Quadrangle arc 2 " << ( aBuilder2.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1367 #endif
1368
1369       VTKViewer_ArcBuilder aBuilder3(P2,P5,P0,myMaxArcAngle); //Build arc using 2, 5 and 0 points
1370 #ifdef __MYDEBUG__
1371       cout << "Quadrangle arc 3 " << ( aBuilder3.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1372 #endif
1373
1374       aCollection.push_back(aBuilder1.GetPoints());
1375       aCollection.push_back(aBuilder2.GetPoints());
1376       aCollection.push_back(aBuilder3.GetPoints());
1377
1378       aScalarCollection.push_back(aBuilder1.GetScalarValues());
1379       aScalarCollection.push_back(aBuilder2.GetScalarValues());
1380       aScalarCollection.push_back(aBuilder3.GetScalarValues());
1381       break;
1382     }
1383     case VTK_QUADRATIC_QUAD:
1384     case VTK_BIQUADRATIC_QUAD:
1385     {
1386       //Get All points from input cell
1387       Pnt P0 = CreatePnt( aCell, inputScalars, 0 );
1388       Pnt P1 = CreatePnt( aCell, inputScalars, 1 );
1389       Pnt P2 = CreatePnt( aCell, inputScalars, 2 );
1390       Pnt P3 = CreatePnt( aCell, inputScalars, 3 );
1391       Pnt P4 = CreatePnt( aCell, inputScalars, 4 );
1392       Pnt P5 = CreatePnt( aCell, inputScalars, 5 );
1393       Pnt P6 = CreatePnt( aCell, inputScalars, 6 );
1394       Pnt P7 = CreatePnt( aCell, inputScalars, 7 );
1395
1396       VTKViewer_ArcBuilder aBuilder1(P0,P4,P1,myMaxArcAngle); //Build arc using 0, 4 and 1 points
1397 #ifdef __MYDEBUG__
1398       cout << "Quadrangle arc 1 " << ( aBuilder1.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1399 #endif
1400
1401       VTKViewer_ArcBuilder aBuilder2(P1,P5,P2,myMaxArcAngle); //Build arc using 1, 5 and 2 points
1402 #ifdef __MYDEBUG__
1403       cout << "Quadrangle arc 2 " << ( aBuilder2.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1404 #endif
1405
1406       VTKViewer_ArcBuilder aBuilder3(P2,P6,P3,myMaxArcAngle); //Build arc using 2, 6 and 3 points
1407 #ifdef __MYDEBUG__
1408       cout << "Quadrangle arc 3 " << ( aBuilder3.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1409 #endif
1410
1411       VTKViewer_ArcBuilder aBuilder4(P3,P7,P0,myMaxArcAngle); //Build arc using 3, 7 and 0 points
1412 #ifdef __MYDEBUG__
1413       cout << "Quadrangle arc 4 " << ( aBuilder4.GetStatus() == VTKViewer_ArcBuilder::Arc_Done ? "" : "NOT " ) << "done !!!" << endl;
1414 #endif
1415
1416       aCollection.push_back(aBuilder1.GetPoints());
1417       aCollection.push_back(aBuilder2.GetPoints());
1418       aCollection.push_back(aBuilder3.GetPoints());
1419       aCollection.push_back(aBuilder4.GetPoints());
1420
1421       aScalarCollection.push_back(aBuilder1.GetScalarValues());
1422       aScalarCollection.push_back(aBuilder2.GetScalarValues());
1423       aScalarCollection.push_back(aBuilder3.GetScalarValues());
1424       aScalarCollection.push_back(aBuilder4.GetScalarValues());
1425       break;
1426     }
1427     case VTK_QUADRATIC_POLYGON:
1428     {
1429       int nbP = aCell->GetNumberOfPoints();
1430       std::vector< Pnt > pVec( nbP + 2 );
1431
1432       for ( int i = 0; i < nbP/2; ++i )
1433       {
1434         pVec[i*2 + 0] = CreatePnt( aCell, inputScalars, i );
1435         pVec[i*2 + 1] = CreatePnt( aCell, inputScalars, i + nbP/2 );
1436       }
1437       pVec[ nbP   ] = pVec[ 0 ];
1438       pVec[ nbP+1 ] = pVec[ 1 ];
1439
1440       for ( int i = 0; i < nbP; i += 2 )
1441       {
1442         VTKViewer_ArcBuilder aBuilder( pVec[i], pVec[i+1], pVec[i+2], myMaxArcAngle );
1443         aCollection.push_back( aBuilder.GetPoints() );
1444         aScalarCollection.push_back( aBuilder.GetScalarValues() );
1445       }
1446       break;
1447     }
1448     default: //Unsupported cell type
1449       return;
1450   }
1451
1452   if(triangulate){
1453     const vtkIdType numFacePts = 3;
1454     vtkIdList *pts = vtkIdList::New();
1455     vtkPoints *coords = vtkPoints::New();
1456     aCellType = VTK_TRIANGLE;
1457     vtkIdType aNewPts[numFacePts];
1458     vtkIdType aTriangleId;
1459
1460     vtkPolygon *aPlg = vtkPolygon::New();
1461     std::map<int, double> aPntId2ScalarValue;
1462     aNbPoints = MergevtkPoints(aCollection, aScalarCollection, aPlg->GetPoints(), aPntId2ScalarValue, aNewPoints);
1463     aPlg->GetPointIds()->SetNumberOfIds(aNbPoints);
1464
1465     for(vtkIdType i = 0; i < aNbPoints;i++) {
1466       aPlg->GetPointIds()->SetId(i, aNewPoints[i]);
1467     }
1468
1469     aPlg->Triangulate(0,pts,coords);
1470
1471     for (vtkIdType i=0; i < pts->GetNumberOfIds(); i+=3) {
1472       aNewPts[0] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i));
1473       aNewPts[1] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+1));
1474       aNewPts[2] = output->GetPoints()->InsertNextPoint(coords->GetPoint(i+2));
1475
1476       if(outputScalars) {
1477         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i)]);
1478         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i+1)]);
1479         outputScalars->InsertNextTuple1(aPntId2ScalarValue[pts->GetId(i+2)]);
1480       }
1481
1482       aTriangleId = output->InsertNextCell(aCellType,numFacePts,aNewPts);
1483
1484       if(myStoreMapping)
1485         InsertId( cellId, aCellType, myVTK2ObjIds, theDimension2VTK2ObjIds );
1486       outputCD->CopyData(cd,cellId,aTriangleId);
1487     }
1488     pts->Delete();
1489     coords->Delete();
1490     aPlg->Delete();
1491   }
1492   else {
1493     std::map<int, double> aPntId2ScalarValue;
1494     aNbPoints = MergevtkPoints(aCollection, aScalarCollection, output->GetPoints(), aPntId2ScalarValue, aNewPoints);
1495     if(outputScalars)
1496       for(vtkIdType i = 0; i < aNbPoints; i++)
1497         outputScalars->InsertNextTuple1(aPntId2ScalarValue[aNewPoints[i]]);
1498     newCellId = output->InsertNextCell(aCellType,aNbPoints,aNewPoints);
1499     outputCD->CopyData(cd,cellId,newCellId);
1500
1501     if(myStoreMapping)
1502       InsertId( cellId, aCellType, myVTK2ObjIds, theDimension2VTK2ObjIds );
1503   }
1504
1505   if (aNewPoints)
1506     delete [] aNewPoints;
1507 }
1508
1509
1510 void VTKViewer_GeometryFilter::SetQuadraticArcMode(bool theFlag)
1511 {
1512   if(myIsBuildArc != theFlag) {
1513     myIsBuildArc = theFlag;
1514     this->Modified();
1515   }
1516 }
1517 bool VTKViewer_GeometryFilter::GetQuadraticArcMode() const
1518 {
1519   return myIsBuildArc;
1520 }
1521
1522 void VTKViewer_GeometryFilter::SetQuadraticArcAngle(double theMaxAngle)
1523 {
1524   if(myMaxArcAngle != theMaxAngle) {
1525     myMaxArcAngle = theMaxAngle;
1526     this->Modified();
1527   }
1528 }
1529
1530 double VTKViewer_GeometryFilter:: GetQuadraticArcAngle() const
1531 {
1532   return myMaxArcAngle;
1533 }
1534
1535 int VTKViewer_GeometryFilter::GetAppendCoincident3D() const {
1536 // VSR 26/10/2012: see description of SHOW_COINCIDING_3D_PAL20314
1537 // in the top of this file
1538 #ifdef SHOW_COINCIDING_3D_PAL20314
1539   return myAppendCoincident3D;
1540 #else
1541   return false;
1542 #endif
1543 }
1544
1545 void VTKViewer_GeometryFilter::SetAppendCoincident3D(int theFlag) {
1546   if(myAppendCoincident3D != theFlag){
1547     myAppendCoincident3D = theFlag;
1548     this->Modified();
1549   }
1550 }